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Abstract 

We formulate a low-energy theory for the magnetic interactions between elec- 
trons in the multi-band Hubbard model under non-equilibrium conditions deter- 
mined by an external time-dependent electric field which simulates laser-induced 
spin dynamics. We derive expressions for dynamical exchange parameters in 
terms of non-equilibrium electronic Green functions and self-energies, which 
can be computed, e.g., with the methods of time-dependent dynamical mean- 
field theory. Moreover, we find that a correct description of the system requires, 
in addition to exchange, a new kind of magnetic interaction, that we name twist 
exchange, which formally resembles Dzyaloshinskii-Moriya coupling, but is not 
due to spin-orbit, and is actually due to an effective three-spin interaction. Our 
theory allows the evaluation of the related time-dependent parameters as well. 

Keywords: Non-equilibrium magnetism; Multi-band Hubbard model; 
Exchange interactions; Non-equilibrium Green functions; Strongly correlated 
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1. Introduction 

Manipulation of the magnetic order of materials via electromagnetic fields 
on increasingly fast time scales is one of the most intriguing issues of modern 
magnetism [1]. While there seems to be a limit (of the order of ~ 2 ps) on 
the time scale attainable in magnetization switching via magnetic fields [2j, pi- 
oneering experimental works dating back to the 1990's demonstrated that 
ultrashort optical laser pulses trigger magnetic phenomena on sub-picosecond 
time scales. These early works opened the field of ultrafast spin dynamics, 
which has been characterized by a flourishing of experimental efforts in the last 
15 years including inverse Faraday effect Q, all-optical helicity-dependent 
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switching of the magnetization |8| , complete write- read events requiring a time 
of only 30 ps [9( , reversal of lattice magnetization in ferrimagnets via a transient 
ferromagnetic state Besides their obvious practical role in the design of a 
new generation of ultra-fast memories, these researches have important impli- 
cations for a fundamental understanding of the interaction between light and 
matter on short time scales. 

From the point of view of theory, the study of magnetization in a realistic 
solid-state system is a challenging problem. In equilibrium, magnetic interac- 
tions in magnetic metals and semiconductors are known to be non-Heisenberg 



1 l!4l4| . that is, the lengths of magnetic moments and values of exchange pa- 



rameters depend on the magnetic configuration for which they are calculated. 
The accurate calculation of these quantities must then follow from an ab-initio 
formulation, and represents a complicated many-body problem. In the case 
of equilibrium, the expressions for computing exchange parameters have been 
given years ago, either within the multiple-scattering formalism in density func- 
tional theory 15[ or, more recently, in terms of electronic Green functions and 
self-energies within the Matsubara scheme, for a multiband Hubbard model 16- 
[l8| . These formulas are commonly used now for a very broad circle of magnetic 
materials, see, e.g., Refs. [l9l - [23 |. 

One could think of using these first-principle methods to calculate equilib- 
rium exchange parameters for further use within a classical Heisenberg model 
to simulate spin dynamics. However, this approach is not expected to allow 
for a satisfactory description of ultrafast magnetism. In fact, the time scale of 
the laser pulses is typically faster than the typical scale of exchange interac- 
tions (« 10 -j- 100 fs). This implies that magnetic interactions cannot be treated 
adiabatically, i.e., the relevant parameters (such as exchange) depend on time. 
Ab-initio spin dynamics [l^,[2g|, also, is based on equilibrium electronic struc- 
ture. 

Within the LDA++ approach 27|, the first-principles electronic structure is 
mapped to the multiband Hubbard model. A multiband Hubbard model, with 
realistic tight-binding and interaction parameters, is likely to be general and 
flexible enough to describe many strongly correlated systems relevant for ultra- 
fast magnetism. In order to include a time-dependent optical excitation, we 
allow the hopping parameters to depend on time. Then, the main approxima- 
tion we can take advantage of consists in the fact that spin dynamics is known 
to be much slower than electron dynamics in relevant systems |3j. This means 
that an effective atomistic model can be derived, with time-dependent param- 
eters accounting for the magnetic interactions mediated by the fast electronic 
dynamics, which can be computed from first principles. Technically, we need to 
separate the spin degrees of freedom from the electronic ones, and derive an ef- 
fective action for the spin variables, after integrating on the electronic variables. 

The needed formalism is that of non-equilibrium Green functions, a very 
general tecnique which resulted from the contributions of several authors (for re- 



views and lectures, see Refs. 28 32]). One typically makes a distinction between 



the equilibrium formalism, due to Matsubara 33[ , and the non-equilibrium for- 
malism developed mainly by Schwinger j34| and Keldysh [3^], which neglects 
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initial (equilibrium) correlations. The Schwinger-Keldysh formulation has been 
applied to the study of spin dynamics, e.g. for a sin g le s p in in a Josephson junc- 
tion [3(| , or in a junction between ferromagnets (371 . l38j , or combined with the 
mean field approximation for the treatment of magnetic interaction [39] . How- 
ever, a first-principle study of an extended and strongly correlated system out 
of equilibrium, such as a fermionic multi-band Hubbard model, has not been at- 
tempted yet, and in this paper we present a general framework for this relevant 
condensed-matter model. Moreover, since we cannot realistically neg lect the 
initial correlations, we need to use the Kadanoff-Baym formalism [28, 3(1 31 1, 
which unifies the approaches of Schwinger, Keldysh and Matsubara. 

The advantage of this approach is that it does not need any assumption on 
the time dependence of the external field, so there is no restriction on time scales, 
which allows to study the role of non-adiabatic and non-Heisenberg effects in 
magnetization dynamics. It is also suitable to make a first-principle formulation 
of quantum noise, whose time scale may be comparable to that of the ultrafast 
pulse, invalidating the Landau-Lifschitz equations, and whose treatment has up 
to now relied on phenomenological parameters [401 ] . 

Besides being an initial step towards a realistic first-principle model of spin 
dynamics in condensed-matter systems, our work is immediately relevant to the 
following problem. Let us consider a magnetically ordered system with negligi- 
ble spin-orbit interactions and anisotropy. Applying a purely electric field does 
not generate spin dynamics, in the sense that it cannot flip the spins of indi- 
vidual electrons, because the electric field does not couple directly with spins, 
and there is no spin-orbit. Nevertheless, a time-dependent electric field mod- 
ulates the possibility of electronic hopping between different lattice sites (and 
orbitals), which in turn affects the strength of magnetic interactions, such as 
exchange. Therefore, if the electron spins are not all parallel in equilibrium, and 
the electronic dynamics is fast enough, then one can assume that the magnetic 
structure, at each time, evolves in the way that minimizes the spin action, with 
instantaneous parameters characterizing the magnetic interactions that depend 
on the time-dependent electronic configuration and hopping. The only con- 
straint on the minimization is that the total numbers of electrons with spin up 
and down along a chosen direction must remain constant. This would be a first 
direct application of the results of this contribution. 

This Article is organized as follows: in Sections 2 - 9 we detail the method 
which we use to derive the effective action for the spins, in Sections 10 - 13 
we present the results, and Section 14 contains our conclusions. In particular, 
Section 2 contains a short review of the Kadanoff-Baym formalism, intended for 
non-expert readers. In Section 3 we present the fermionic multiband Hubbard 
model, and in Section 4 we show how to transform from the fermionic represen- 
tation to a bosonic one, where the bosons are related to the directions of the 
spin axes. In Section 5 we discuss how to integrate away the fermionic degrees of 
freedom from the total action, obtaining an effective bosonic action for the low- 
energy excitations on the top of the ground state, which is shown in Section 6. 
In Section 7 we introduce the Kadanoff-Baym equations of motion, which allow 
to express the magnetic interactions appearing in the bosonic action in terms 
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of fermionic Green functions and self-energies, as we do in Section 8; in Section 
9 we show a way to slightly simplify the action. After this introductory part, 
in Section 10 we start to present the results of the derivation, showing that the 
out-of-equilibrium multiband Hubbard model, close to the equilibrium ground 
state, can be mapped into a rather complicated effective model that includes 
two kinds of time-dependent magnetic interactions. In Section 11 we show that 
our formalism correctly recovers the known results [lfj in the equilibrium case, 
while in Section 12 we discuss the general non-equilibrium case, together with 
some relevant approximations which can simplify the calculations. In Section 
13 we give the expression of the dynamical spin stiffness, which is an interesting 
observable in view of comparison with experiments. In Section 14 we summa- 
rize our results and discuss possible further developments. Appendices A - G 
include mathematical details. 



2. Kadanoff-Baym approach to non-equilibrium systems 



For the sake of completeness, we here review briefly the main concepts of 
the Kadanoff-Baym approach to the study of non-equilibrium systems 28l-30j. 
which combines the approaches of Matsubara [33| , Schwinger [34| and Keldysh 
23- 



2.1. Kadanoff-Baym time contour 

Given a time-dependent Hamiltonian H(t), the equation of motion for the 

density operator is dp(t)/dt = — i H(t),p(f) , which can be solved formally 

as p(t) — U(t,to)p(to)U(to,t), where p(to) is the (supposedly known) density 
operator at a reference time to, and the evolution operator is 



U(t,t') = e{t-t')Texp 



1 J dhHih] 



+ 6(t' - t)Texp ( -i I dtiH(h) 



(1) 



for t 7^ t 1 , while U(t,t) ~ 1; the symbol Q(t) denotes the step function [Q(x > 
0) = 1, Q(x < 0) = 0]. The expectation value for the observable O at time t is: 



Tr 



0(t) 



6p(t) 



Tr 



p(t )U(ta,t)dU(t,t ) 



(2) 



Tr [p(t)\ Tr [p(t )} 

where the trace is evaluated over the complete many-body Hilbert space, and in 
the last passage the cyclic property of the trace and the identity U(t, t')U (t', t) = 
U(t,t) — 1 have been used. We choose the reference time to in such a way 
that, for t < t , the Hamiltonian is independent of time and the system is in 
equilibrium. Therefore, we can use for p(to) the grand-canonical equilibrium 
expression, 

-f)(H -fi,N) 

p(*o) = r o/ -, — , (3) 



Tr e-PiVo-vN) 
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where Hq = H(t < to), N is the numbcr-of-particle operator and /j, is the 



chemical potential. We assume 



H ,N 



= 0. We extend the time domain to 



the complex plane, defining the complex time variable £ = t — it, with the 
understanding that H (() = H(t) depends only on the real part of time, and we 
define the evolution operator in imaginary time 



e-W°-»*)=U v (t -U3,t )- 
Thus, we can write Eq.([2]) as: 



(4) 



Tr 



0(t) 



U v (t - i/3, t )U{t , oo)U(oo, t)OU (t, t Q ) 



Tr 



a -/9(j& -|*#) 



Tr 



U v {t Q - i/3, t )U(t ,t)OU(t, oo)*7(oo, t ) 



Tr 



-0(H o - f j,N) 



(5) 



which shows that the evaluation of 0(t) requires that we let the system evolve 
along three time domains in the complex plane: a forward branch 7 + = (t , oo), 
a backward branch 7_ = (cxD,t n ), and a segment on the imaginary (vertical) 
axis of time 7„ = (to, to — i/3), which we call the vertical branch (see FigflJ. It 
must be noted that the time value which we have labelled as oo may actually 
be chosen as a completely arbitrary (finite) value. The total domain over which 
the system evolves is the Kadanoff-Baym contour 



r = 7+ U 7 -U7„ = (io,oo)U(oo,i )U(i ,t -i/3). 
We define the total evolution operator on T as: 

U T = U v (t -i/3,t )U(to,oo)U(oo,t ). 



(G) 



(7) 



From Eq.JS]), we see that the computation of 0(t) is realized by opening the 
contour T at the instant t either on the branch 7+ or on the branch 7^ , inserting 
there the Schrodinger-represented operator O in Ur, and evaluating the trace 
of the resulting operator. The inclusion of the branch 7^ is required to treat 
systems where the inital correlations are not negligible, as it is the case in 
typical solid-state systems. In the cases where the initial correlations can be 
neglected, the initial density matrix is of the single-particle kind and there 
is no need to express it in a contour formulation: in such conditions, one can 
restrict the contour to 7+ U7_, possibly with to — > —00, which is the Schwinger- 
Keldysh contour. In equilibrium conditions, on the other hand, the Hamiltonian 
is time-independent and the contour is restricted to 7l) , which is the Matsubara 
contour. The Kadanoff-Baym formulation, therefore, unifies and generalizes 
the other approaches, and allows to treat the most general case of a system in 
non-equilibrium with initial correlations. 
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Figure 1: Schematic representation of the Kadanoff-Baym contour. Branches 7-1- and f— are 
displaced for graphical convenience, but they both coincide with a portion of the real time 
axis, extending from to to 00. 



2.2. Path integral formulation 

The Kadanoff-Baym partition function is defined as 



Tr 



'Tr 



-0(Ho-hN) 



(8) 



where Ur is given in Eq.©. Since U(to, oo)C/(oo, to) = 1, it follows that Z — 
1. Despite this apparent triviality, expressing Z by means of path integrals 
allows to derive a non-equilibrium action, from which it is possible to extract 
physical information (32)]. To do this, we start by denoting the Hamiltonian as 

H (f>',(j);t for t > to, where $ and (f> represent the sets of fermionic creation 
and annihilation operators, respectively. Instead, for t < to and on the branch 



7„ the Hamiltonian is constant and we denote it as Hq 
the branch j v by means of the real variable 

r = -Im(C), 



We parameterize 



(9) 



which is equal to for £ — to and to j3 for C, = to — i/3. On the 7^ contour, 
£ = to + i I m (C) = ^0 — i r - Then, standard manipulations lead to the expression 



V [<f>, 0] e iS [^] 



(10) 



where the effective action S [0, (f\ , written in terms of the Grassmann variables 

operators, is: 



relative to the 

S\$,6\ : 



*0+e 



+ i 



t{i4>+(t) ■ Mt -e)-H [<£+(*), <M* - e);t] 
i^_(i-£)^.(t-e) + ff [0_ (t-e), 4>- (t); t] } 
dr{^(r) -4> v (T-e)+K [&(t),&,(t - e)] }, (11) 
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where K = Ho — fiN, and which requires some explanations. First, Eq. flTTT) 
is written, for convenience, in terms of real time variables t and t, instead of 
contour variables. Therefore, since each value of t corresponds to two distinct 
points on the Kadanoff-Baym contour (one on 7+ and one on 7-), the time- 
dependent Grassmann fields must be specified by the index ± if their argument 
is on the real-time axis, while we use the label v for the fields with the argument 
on the branch j v . Then, we have introduced two infinitesimally small positive 
quantities, e and e, in order to emphasize a subtle technical point, namely the 
fact that any product of $ and <fi fermionic operators appearing in the Hamil- 
tonian transforms in the path integral formulation into a product of <j) and 4> 
Grassmann fields in which the fields <f> are evaluated at an instant occurring 
infinitesimally before (in the contour sense) the instant when the </> fields are 
evaluated [32J]. While this aspect is often neglected, on the basis that the fields 
are assumed to be continuous functions, for our purposes it will be important to 
keep explicitly track of this discrete structure because we will have to deal with 
discontinuous functions, such as the correlators (Green functions) originating 
from the Grassmann numbers ((/>,</)). In these cases, we will need to consider 
carefully the direction along which the independent variables approach the dis- 
continuity point, which will be possible in our formulation. At the appropriate 
stage of the calculations, we will send e — > + and e — > + . On the other hand, 
whenever e or e appears inside the argument of a continuous function, we are 
allowed to send it to immediately, as we already did implicitly in the case of 
the time-dependent external field. Finally, derivatives of fields are always meant 
to be taken from the right side, i.e., <j>(t) — lim e _j. + [<j){t + e) — 4>(t)] /e. 



3. Multi-band Hubbard model 

In order to model an electronic system driven out of equilibrium by a time- 
dependent external field (e.g. a laser pulse), we must consider a Hamiltonian of 
the form 

H(t) = H T {t)+H v , (12) 

where Hxit) is the single-particle Hamiltonian, including the time-dependent 
field, and Hy is the (time-independent) interaction potential between the elec- 



trons. We will treat a multi-band Hubbard model [27|, [41] , therefore the elec- 
tronic single-particle states are identified by three labels: the site index i, the 
orbital index A and the spin index a. We ignore spin-orbit coupling and assume 
that the external field is diagonal in spin indices (we are therefore excluding 
magnetic fields, but including purely electric fields which are relevant for mod- 
elling all-optical experiments). The single-particle Hamiltonian is then given 

by 

HHt) = E r <.A.wt) E ^La4v = E E T -"(t) E 

i a A n £(,A& c a b <y 

a b 
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where we have grouped the site and orbital indexes according to a = (i , A a ) 
and 6 = (ife, Aj,), and we have defined the spinor fermionic operators 

«-(«, (") 

The matrix element T ab {t) — T b * a (t) of the single-particle Hamiltonian is written 

as 

T ab (t) = T ab + f ab (t), (15) 

where T ab is the time-independent hopping parameter due to electronic struc- 
ture, and fab(t) is the time-dependent matrix element of the perturbing field. 
We denote as to the time at which the external field is switched on: f a b{t) = 
for t < tg. The interaction potential generating Hy is assumed to be on-site, 
i.e., 

t A1A2A3A4 era' 

Including magnetic fields and spin-orbit coupling complicates the analysis 
significantly, and will be left to future work. However, this model already al- 
lows to describe some interesting magnetic phenomena. For example, consider a 
system which contains both spin-up and spin-down electrons. The arrangement 
of the spins within each lattice site depends initially on the equilibrium mag- 
netic interactions, mainly exchange. Applying a time dependent electrostatic 
field on a portion of the sample, as we shall show, may change the strength 
of the magnetic interactions. If in a certain region of the sample the coupling 
switches, e.g., from antiferromagnetic to ferromagnetic, then this might gener- 
ate a re-arrangement of the total spin in each lattice site as a 'purely electronic 
phenomenon, due to the inter-site hopping. Hence, even without individual spin 
rotations, domains with an ordering which is different from the initial one may 
originate as a consequence of electronic transfer between atomic sites. 



4. Rotation of the spin quantization axes 

4-1. Holstein-Primakoff bosons 

For each site, we define a rotation matrix, acting in the space of spinor 
fermionic operators, as: 



r IMV " , * w j > (17) 




where z is a contour variable which parameterizes the Kadanoff-Baym contour. 
In Eq. pTl) we have introduced the boson fields 

^^-e 1 ^ sin [0^/2], (18) 
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with Oi 6 [0, 7r[ , y>i € [0, 27r[ being the polar angles that determine the spin axis 
on site i at time z; it holds that R\(z) ■ Ri(z) = 1. 

We transform the Grassmann variables appearing in the action, Eq. pT]) , 
according to 

4>a±{t) = $ a± (t) ■ Rl ± (t), 4> a ±{t) = R a± (t) ■ l/> a ±(t), 

4>av(r) = ^av(r) ■ RI v (t), (f> a v(r) = R a v(r) ■ ipavir) ■ (19) 

To understand the meaning of the rotation that we have just introduced, the 
following considerations are in order. The local vector spin operator for site- 
orbital a, in the laboratory reference frame, is 

& a = ft a -a4 a , (20) 

where <x is the vector of Pauli matrices. The expectation value of this operator 
on the state d>t„\0) is 



<JU Z , 



(21) 



where a € {ti4-} = {+? ~ }■ Instead, the expectation value of the spin operator 
on the state ^y aa |0) is 



(22) 



where we have used Eqs. (fT51) and <jTUJ) , and the unit vector e a is given by 

e a = u x sin(0 o ) cos((p a ) + u v sin(6» a ) sin(^ a ) + u z cos(6» a ). (23) 

Therefore, e a (which depends on time) has the meaning of the unit spin vector 
on site-orbital a, measured in the laboratory reference frame, if the site a hosts 
a -0-)- electron 42] . Expression (|2"3"|) can also be written as 



(24) 



which shows that ^ and ^* are, effectively, nothing else than Holstein-Primakoff 
bosons SJfil. 



4-2. Transformed action 

We assume the rotation matrices R(t) to be differentiable functions of t 
over the interval ]tg + e, oo[, and the rotation matrices R(t) to be differentiable 
functions of r over the interval ]e, (3[. So, in the evaluation of Eq. pTj) we can 
already apply the replacements R(t — e) — > R(t), R(t — e) —> R(t). Instead, we 
need to keep finite e and e in the arguments of the fermionic fields (i/j,ip). 

Since we are considering a on-site interaction, Hy is rotationally invariant, 
hence it is unaffacted by the transformation defined in Eq. (|19[) . which means 
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that Hv[4>,<p] = Hv[ip,ip}- On the other hand, the single-particle Hamiltonian 
acquires a dependence on the bosonic fields, Ht[4>, <p] t] = Ht[$, tp,£* ,£;t}. 

The transformed action, depending on the fermionic (Grassmann) ip fields 
and on the bosonic (complex) £ fields, is written as: 



s [0, 4>]=s [i>, i>,?,t]=s + s' v, r, C] 



(25) 



with S Up, ip\ denoting the original expression of the action, with the (<jf>, </>) 
fermions replaced by the fermions, and 



s' [$,ii>,e,t] 



dTip v (r) ■ A„(t) -ip v (T- e) 



/•DO 

- / dt 

Jto+e 



4>+(t) ■ A+(t) ■ i/>+(t - e) - - e) • A_(t) • ^_(t) 



(26) 



where the quantities A„ and A± are matrices in both the spaces of Hubbard 
indexes and spin indexes, and they depend on the fields and (£±,£±), 

respectively. Expressing explicitly their Hubbard-space structure, 



A„6,±(t) = iR f a± (t)R a± (t)6 ab - T ab (t) [Rl ± (t)R b ±(t) - 1 
( r ) = iRl v (T)R av (T)S ah + \T ab [rI v (t)Ri, v (t) - l] . 



(27) 



These are now matrices in spin space, which include diagonal and non-diagonal 
components. We separate the corresponding parts of the action: 



with 



EEE 

a b <y 

/'OC 

+ / dt 

Jto+e 

EEE 

a b <y 

/'OC 

+ / dt 

Jto+e 



S' = + £(2) 



(28) 



dC(T)A:yr)C(r- £ ) 
A£ + (i) - 6) A£_(t) V£_(t) 



(29) 



where we have introduced the notation a = —a. 

4-3. Partition function 

The partition function is written in terms of the rotated fields as: 



(30) 
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where the symbol 



(31) 



means integration over all possible orientations (8 a , fa) at all lattice sites. The 
dependence of the fields (£a,£ a ) on the angles (9 a ,f a ) is given by Eq.([Tg)). It 
can be shown that integrating over the angles (9 a , f a ) ensures that the partition 
function as written in Eq. (l30l) is equal to Eq. lfTC))) . i.e., Eq. (f30|) is an identity. 
Since f is the phase of the complex number — £, while sin(0/2) is its modulus, 
we see that 



271 



d0sin(0) / dff -sin(6>/2)e IV ,-sin(6l/2)e- IV 



= *jj o dRe(0 dim® /[£,£*], 



(32) 



where the integration domain C\ is the circle of radius equal to 1 in the com- 
plex plane, centered on 0, i.e., |£| 2 < 1, which is described by f <G [0, 27r[ and 
9 € [0, 7r[. We therefore change the path variables from (9,f) to (£,£*), by 
introducing the notation 



so that 



Jv[d,f\=l[ ^JJ c dRe(Ca) dlm(Ca) = J 



(33) 



(34) 



4-4- Low-energy theory 

The introduction of the Holstein-Primakoff bosons allows to decouple the 
dynamics of electronic and spin degrees of freedom. Now suppose that the 
equilibrium spin configuration is collinear (ferro-, antiferro-, or ferri- magnetic), 
with u z being the initial direction of atomic spin alignment. Our goal is to 
derive a theory for the low-energy excitations on the top of the equilibrium 
configuration. Such excitations corrispond to small deviations of the atomic 
spins from the direction of u z , i.e., we can assume that the polar angles 9i be 
small. We can therefore approximate 



1 



5« 1 2 / 



(35) 



Now the unitarity of the transformation (|19|) holds but for corrections of the 
order of |£i| 4 . Under this approximation, 



7?t R 



Sas 



(36) 
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R\ s Rbs — 1 ~ f 



£as ^bs 



^bs ^>as 

!€a £ l!+|6, 



(37) 



This procedure is justified if we limit ourselves to the description of small rota- 
tions of the spins with respect to the direction of the initial quantization axis. 

For future convenience, we define permutation operators P ab , which substi- 
tute a with b and b with a in the expressions they act upon, where a and 6 
are Hubbard indexes. Moreover, in order to adopt a compact notation, we will 
sometimes put arrows — >, <— over the permutation and derivative operators to 
indicate the direction along which they act. The A matrices are written as: 



AJL(*) = 

Alt(r) 
AtL(r) = - 



S ab ig- t +T ab (t) (l-f al 
■it 



C±(*)> 

6»,±(*)> 



A^ ± (t) = ^a±w^ i ^-^ )C±W 



abi ^ - iT a6 (l - ~P ab ) CM, 



+ T o6 (t) 



Ait i± (t)=J„6^ ± (^(^-^ka ± (t) 



2\dt dt 

\Za±(t)\ 2 + \t b ±(t)\ 2 

2 



+ T a6 (i) 



2 \ dt dt 

\Za±(t)\ 2 + \(; b ±(t)\ 2 



AH, t (r)=«aA)5|i-i]ar) 



<9 



iT, 



a6 



Att(r)=^CW^[i--^-)U(r) 



2 I <9r <9t 

\L v (r)\ 2 + \^ bv (r)\' 



i ( d d 



-iT 



0.6 



2 \^dr 0r 

|U(r)| 2 + |^(r)| 5 



- Lv(r)H* bv (T) 



(38) 



(39) 



(40) 



(41) 
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It is seen that the Aj^'s are linear, while the A^'s are quadratic, in the bosonic 
fields (£*,£)• 

To be consistent with the small #'s approximation, we must expand the 
partition function in series of the bosonic fields and retain only up to quadratic 
terms. We obtain 



Z = j V $,ip,e,(\ e i { s I'M+ s '[^>«*>«]} 

« Jv[^,C,(\ e iS ^{l + iS' -\{S^ $,tl>,C,£] 

= Z<®+ZV>+ZV\ (42) 
where 

Z^ = [ V i>, e iS Mi<?( 2 ) ^ (43) 



5. Fermionic correlators 

In the next sections, we will integrate out the fermionic fields, by exploiting 
the expression for the single-particle Green function (for s 6 {+,—,?;}, and 
z = t or z = t, accordingly), 

V [j,,1>] e iS &*\$£,(z')il&(z) ee G%Z' ss ,(z,z'). (44) 

Since the Hamiltonian is spin- independent (i.e., the hopping is diagonal in spin 
space), the Green functions with a ^ a' are zero in our system. Therefore, 

V[^,CA e iS [<M i5 (i) =0 

=> = |l> ,£] e iS M {-- (SM $,il>,?,t]) 2 } ■ (45) 

As a consequence, the action will contain no terms linear in the bosonic fields, 
which is due to the fact that we are not including spin-orbit coupling. Therefore, 
our model will not include Dzyaloshinskii-Moriya interactions [18| . Since we will 
need only spin-diagonal single-particle Green functions, we will label them with 
just one spin index. We also need the expression for the following two-particle 
Green function, 

(46) 
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We adopt the approximation 



Xbab' a' ,sss' s' i z > z i z i z ) 



-G 



ba' ,ss 



(47) 



which corresponds to neglecting the vertex in the two-particle Dyson equation, 
and where we have already taken into account the fact that correlators G' 7 ' 7 are 
zero in our system. It must be noted that Eq. (|47p is the only approximation 
that we adopt on the many- body level. For the equilibrium case, the same 
approximation was used in Ref. |17j and it was shown that it leads to the standard 
expression for exchange parameters 15, 16| . 

All quantities will then be written in terms of single-particle Green func- 
tions. The correlators can be classified according to the positions of their time 
arguments on the Kadanoff-Baym contour, and put in correspondence with non- 
equilibrium Green functions written in terms of the field operators, as follows: 



G. 



ba,++ 



(ti,* 3 ) = G&r(ti,t a ) 



i(T t 



Gl a , + .{tx,h) = G<:{t l ,h) = y{^l a {hW b {ti)) 

G a ba ,-+(h,t 2 ) ee (?>'(*!, t 2 ) = -i(^(*i)^'(* 2 ); 
Gba, vv (n,r 2 ) ee Gf:{ n - r 2 ) = -i (% ^(n^'fa) 
Gl a>±v {t,r) = G&(t,r) = i , 
G° ba<v± (T,t) = G r C( T ,t) = -i(#(r)#"(i)) , 
where we have used the notation 

(d(z u z 2 ))=Tr{e-f , ( 6 °->**)d(z 1 ,z 2 )} /lY {e^o-^)} , 

and 

4>(t) = U{t ,t)TpU(t,t ), 

V3(r) = U v (t Q ,t - iT)ijjU v (t - ir,t ) = A^-^)^e-< 6 "-^), 



(48) 



(49) 



(50) 



where the identity U v (to,to — it) = c T ( J ^ 0_A1 ^ r ) follows from the fact that the 
Hamiltonian is constant (equal to Hq) on the vertical branch of the Kadanoff- 
Baym contour. The correlators G^{t\,t 2 ) and G^(ti,t 2 ) are continuous func- 
tions of their time arguments (ti,t 2 ), while Gj^(ti,t 2 ) and Gj%(ti,t 2 ) are dis- 
continuous across the line t\ = t 2 if b = a. In fact, from the definitions of the 
Green functions related to the real-time branches, it follows that, for t 7^ t' , 



gL(m') = e(t - t')G> a (t,f) + e(f - t)G< a (t,f), 

Gl(t,t>) = &(t> - t)G> a (t,t>) + Q(t - t')G< a {t,t'), 



(51) 
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where Q(t — t') is the step function. Therefore, 

lim Gj a (t, t + e)= lira G T ba (t - e,t) = G< (t, t), 

c^0+ e— S-0+ 

lim GL (t, t - e) = lim Gj, (t + e, t) = G> (t, t), (52) 

and 

lim Gl(t, t + e)= lim Gf Q (f - e, t) = G> (t, i), 

£— >0+ e— S-0+ 

lim GL (t, t - e) = lim cf Q (t + e, t) = G< (t, i). (53) 

c— >0+ e— S-0+ 

The discontinuity is: 

G<(M)-G>(M) = i<U. (54) 

We introduce the notation 

Pba (t) = -iG<(t,t) = = (u(t ,t)i>li> b U(t,t )) (55) 

to denote the time-dependent pair correlation function; in the case a = b we will 
drop one of the identical subscripts and denote as p a (t) the occupation number 
on site-orbital a at time t. From the definitions, Egs.(|48p. it is seen that the 
following property holds: 

'Gf b (t,t')]* = -Gf a (t',t), (56) 

Since the Hamiltonian is constant on the vertical branch of the Kadanoff- 
Baym contour, the Matsubara Green function depends on (t, t') only via the 
difference (r — r'), as indicated in Eqs. (|4"5]) . Therefore, the Matsubara Green 
function G M (r) depends on one time argument only, r € (—f3,f3). It is conve- 
nient to put 

G2(t - t') = 0(r - t')GZ (t - r') + 6(r' - r)G« (r - r'), (57) 



which is consistent with the definition of G given in Eqs. (|48)) . One can imme- 
diately see that 

lim GZ(-e) = Gg,(0) = i U^ b ) = G<(t ,t ) = \ Pba , 

ton Gf a {e) = GZ (0) = -i = G> a (t , t ) = i( Pba - S ab ), (58) 

where we have introduced the equilibrium pair correlation function, i.e., the pair 
correlation function at t = to, p ba = — iG^(to, to). Therefore, 

lim [G&(-e)-G&(e)] =i6 ab . (59) 

e— >0+ 
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Values of G m (t) at the borders ±/3 and close to are related: 

G M (/3) = — G M (— e), G u {~/3) = — G M (e), (60) 

where e — >• + . 

Finally, the functions G n and G r are continuous, and they satisfy the fol- 
lowing boundary conditions: 

Gbafo, 0) = i (i>ti>b) = iPba, 

G r ba (0,t a ) = -i(^l) - i( Pba - S ab ). (61) 

In the following, we will express the partition function in terms of an effective 
action involving the (£, £*) fields only, which will be achieved by integrating the 
fermionic fields out, by means of Eqs. (|48|) . We will apply the Keldysh rotation 
to the bosonic fields living on the Schwinger-Keldysh contour [32| . 



&+(t)\ J_(l 1 \ (£ac(t) 



(62) 



5.1. Evaluation of Z^ 

We now calculate Z^ explicitly. From Eq. (|4"51) . we see that this quantity 
requires the evaluation of double time integrals only. Since the domain where the 
Green functions G(z,z') is discontinuous is one-dimensional (the line z = z'), 
so it has zero measure with respect to the integration domain, we can already 
take the limits e — > + and e — > + . This gives: 

J a, a' 6,6' I Jt ° Jt ° 

( Gg,W)GZW,t) -G<^(t,t')G>l(t',t)\ (A$ Vt+ pj\ 
\-G> a \(t,t')G<l(t',t) G^,(t,f)G^ a (t',t) J \&# V ,-V)J 

dr f dr'Alt (r) A^.Jr') G#( T ~ r>) G^(r' - r) 





oo 

" ' A%Jr)G b l,{t,T)G b 7 a {T,t) 



f°° C@ r 

/ dt dT^A^(i)-At(t) 

J to JO „ L 



(63) 



We now perform the Keldysh rotation. Since the A± 's are linear functionals 
of the £±'s, they transform under the Keldysh rotation exactly as the fields 
themselves, i.e., 



(A& + (t), Al{_(tj) = (Al{ c (t), A%(t) 



M 1 -v 



(64) 
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Therefore, the first term in Eg. (155)) produces a coupling between the £c and the 
£q fields via the matrix 



11 1 

2 \1 -1 







G% (t, t')G T b } a (t\ t) -G<i(t, i')G> T a (i', t) N 



1 1 

1 -1 



(65) 



whose non-zero elements are: 

X* 6a , (t, == G&(t, G>I(i', t) + G>i(t, f) G<£(f , t), 
X§ aba ,(t,t') = 0(f t) [G<j,(t,t') G>l(t',t) - G>i(t,f ) G#[(f ,t)" , 

x « o6o ,(t,f) = -e(t_t') ^(t.OG^il-G^O^C.i)] ■ (66) 

The third term in Eq. (|63| produces a coupling between the £q and the £ v fields, 
since A£f >+ (t) - A££_(i) = V2A£f )Q (i). Thus, gives no coupling at all 
between the £c and the £ v fields. Combining the three terms, the resulting 
expression is 

= - / I'M EE { /°° d< f ^ ( At icW> A&gM 

J a,a' 6,6' I Jt ° *° 

' \x§ aba ,(t,t>) x $ aba >(t,t')J {&$ VtQ V)) 



dT 



dr'AHjr) G#(t - r') A^ b , >') G^(r' - r) 

r/3 



+ v / 2 / di / dr E ^g(t)C(^)A^,„(r)G 6 r f fl (M) 



(67) 



5.5. Evaluation of 

We now calculate Z( 2 ) explicitly: 



J a b a { J z 



dt 



GJZ(t - e , t) A^,+ (*) - G^(t, t - £) A^_ (t) 



(68) 



where we have used the identity G^ vv (t - e, t) eee G^(-e). From Ea.(|57). it 
follows that 



Hm + CN = GSr(0). 



(69) 
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From Eqs. (f52")l and (f53"f we see that 

lim Gl: (t - e,t) = lim Gf:(t, t-e)=i (u(t ,t)^^U(t, to)) = G< ff (t, t). 

(70) 

These substitutions remove the discontinuous correlators and we are thus al- 
lowed to put e = and e = in the integration limits in Eq. (|68l) . This gives: 



J a b * l Jt ° 



dtG<:(t,t) [A£ + (t)-A£_(i)] 



-OO) [" drA% Jt)\. (71) 



We now perform the Keldysh rotation, defined in Eq. (p2^1 . The A c " 7, s are not 
linear in the £ fields, differently from the A c " 7, s, so they do not transform linearly. 
The result for Z^ is: 



J a b a [ Jt < 

• £ 



to 



T, 



^5(ui'.)i,i„i') ■ u(')Ci') 



'« t (^u(c s w - T ab (t) (5 a ^ ab + ^) umut) 

+ G%(0) [" dr &%Jt)\, (72) 



where we have defined the Keldysh index n € {C, Q}, and with k we denote the 
complementary of k. 

6. Effective action 

The three parts of the total partition function can be written, respectively, 

as 

(73) 

where S\ and £2 are effective actions which depend functionally on quadratic 
combinations of the fields (£*,?)■ We can put 

Z = J 2? (l + +iSa J 2?[r,e]e i5r ' cl , (74) 
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where the effective action is S [£* , £] = Si [£*, £] +52 [£* , £] . This total action can 
be written as the sum of several terms, all quadratic in the Holstein-Primakoff 
fields, which we classify according to which fields they couple. We label as S xy 
the action term coupling the fields £* and £ y , with (x, y) € {C, Q, v}. Si provides 
terms Sqq, Scq, Sqc, S vv , Sq v , S v q, while S 2 provides terms Scq,Sqc,S vv . 
The total action is written as: 



S — Sqq + Scq + Sqc + S vv + Sq v + S v c 



(75) 



The individual terms, including contributions from both Si and St, can be 
presented as follows: 



iS QQ 
iScQ 
!Sqc 

i$vv - 
iSq v -- 
iS V Q -- 



dt 



dt' &(t) ■ A(t,t') ■ 



it, 



dt\e Q (t)-C^{t)-ic{t)+ / At'C Q {t)-C^\t,t')-ic{t') 



to K J to 

^ dr [C(r)-D {1) {r)-Ur) + l" & ■ D^(t,t') ■ , 







dt / dTC Q (t)-E(t,T)-UT) 



It, 



dt / drC(r)-F(T,t)^ Q (t), 



(76) 



where all quantities are either vectors or matrices in Hubbard space, and we have 
introduced several kernel operators coupling the bosonic fields, which contain 
time derivatives acting either on the fields at their right side or at their left side. 
Such operators can be expressed conveniently by using the property that the 
arrow over a permutation operator can be flipped under the sign of summation 
over both its indexes, i.e., 

y fabP ' ab9ab = ^ ^ fba9ab = ^ ^ fabfjba = ^ ^ fabP abdabl (77) 



a b 



a b 



a b 



a b 



this property allows to write the individual action terms in Eq. (|76p as matrix 
products in Hubbard space. To simplify the notation, for each spin-dependent 
function f a we define the charge combination, f c = (ft + f^j /2, as well as the 
spin combination, / s = (/^ — /^j /2. We further use Eqs. ([55l) and (|55|) . The 
expressions of the kernel operators are then given by: 



6 b' 



Sabi^+(l-^ab)T,t h (/) 



Xb'aba' (*>*') 



5a'b'i^+T a > V {t') (l 



(78) 
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1 



B^it) =6 aa ,^-p s a (t) - P l(t)-+iJ2 [p° b (t)T ba (t)+T ab (t)pZ(t) 

pi>a(t)Taa'(t)+ P i a ,(t)T a , a (t)} , 



(79) 



b V 



<WtT7 + 



dt 



(l - ~P ab ) T ab (t) 



X b >aba'{t->t') 



(80) 



(81) 



b 6' L 



4- 



abi di ^ ab ) 



XbLba'(M') 



(82) 



(r aa , P \, a + P l aa , r a , a ) , 



(83) 



b V 



abig^ - (l - ^ab) lT ab 



E aa ,(t,T) = V2j2J2 



b b< 



8 ab i^- t + (l - ^ab) T ab (t) 



(84) 



Gli(t,r)G b J a (r,t) 



(85) 



F a , Q (r,t) = ^EE 



b b' 



<WiJ^ - (l - ^a'b') iT^'b' 

abi^+T a6 (i) (l-Kb) 



G^(t,r)G^,(r,t) 



(86) 
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In all of the above expressions, the action of permutation operators extends out 
of the square brackets (so, one first has to eliminate the square brackets by doing 
the products, and then applies the permutation operators in the appropriate 
directions). 



7. Kadanoff-Baym equations 



7.1. Green functions and self- energies 

We now use the Kadanoff-Baym (KB) equations of motion (EsOM) for the 
Green functions [3l| to simplify the expression of the action. We need the 
following relations and definitions: 



X H -G> 



G>{t,t') ■ G"V) 



G^ ■ Y, A + G R ■ 



£*G 



(M') 



~G*-\t)-G\t,T) = [£ K -G n + £ n *G M ] (t,T) = Y\t,r), 

G r (r, t) • G"- 1 ^) = [ G r • S A + G M * E r ] (r, t) = Y r (r, t), 

G\t,T)- t^(r) = [G R • + G" * E M ] (t, r) = r), 

V-^r) • G r (r, t) = [E r • G A + S M * G r ] (r, i) = Z r (r, i), 

f- 1 ( T )-G M (T - t') = \8{t - t') + [S M * G M ] (r, r') = i<5(r - r') + 7 M (r, r'), 

G M (r - /)• V~V) = W( T - r') + [G M * E M ] (r, r') = i<5(r - r') + J M (r, r'), 

(87) 

where all the quantities are matrices in the Hubbard space, and the matrix 
elements of the operators ~ 1 (t), G _1 (i), ~~$~ 1 (t) and V _1 ( r ) are: 



3 

-S ab i— - T ab (t), 



9t 



^ab — (Tab — M<5ab), 



OT 



(88) 



In Eas. (l87l) the non-equilibrium self-energies (S > , £<, £ , £ r ) and the equilib- 
rium (Matsubara) self-energy £ M have been introduced, as well as the retarded 
(R) and advanced (A) combinations of Green functions and self-energies: 

G R (t,t') = 0(t-t') [G>(t,t') - G < (t,t')] , 

G A (t,t') = -6(t' - f) [G>{t,t') - G<(M')] , 

£ fl (t, t') = S(t - *')£(*) +e(t-t') [£> (t, t') - £< (t, 0] , 

£ A (i, t') = 5(t- t')E(t) - 6(i' - t) [E> (i, i') - £< (i, t 7 )] , 
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where E(t) is the Hartree-Fock self-energy 3l|. We have also used the following 
notations: 



/■oo 



[f-g](z,z') = / dtf(z,t)-g(t,z'), 



[f*g](z,z') = -i / drf(z,r)-g(T,z') 



(90) 



o 



From Eqs. (|87j) and Eq. ([56|) . it can be proved that 



(91) 



7.2. Symmetry breaking 

Despite the fact that the Hamiltonian is invariant under uniform rotation of 
all spins around a same arbitrary axis, the ground state exhibits a spontaneous 
symmetry breaking (SSB). The aim of our treatment is to treat low-energy ex- 
citations on the top of the magnetically ordered ground state. The direction of 
the preferred quantization axis of the ground state, u z , may be selected by mi- 
croscopic magnetic fields, which are not included explicitly in our Hamiltonian. 
To implement the SSB, we assume that Green functions are spin-dependent [for 
simplicity, the spin index was not written explicitly in Eqs. (|87l) ]. with G^ ^ G^, 
where the source of the spin dependence is included in the self-energies, E^ 7^ E-k 
In reality, this corresponds to spin-polarized self-consistent solutions of density 
functional [l5| or dynamical mean-field theory [l6| equations. 

7. 3. Properties of the Matsubara functions 

From Eqs. (|87|) . using Eqs. (|57l) and (|58|) . one obtains the equations for the 
continuous functions G u and G D , valid for t^t': 



0(t - r') f-\r) ■ G a (r - r') + Q(t> t) ^(t) • G v (r r') = 7 M (r, r'), 
G D (r - r') • V~ V) 6(r - r') + G u (r - /) ■ V~ V) ©(^ - r) = J M (r, /)■ 



Moreover, from the definition of the Matsubara Green function, the following 
boundary conditions follow: 



(92) 



G M (r + /3) 
G M (r - p) 




for - /3 < t < 0, 
for < r < p, 



(93) 



which can be written compactly as 



G M [r -sign(r)/3; 




for - (3 < t < P A t ^ 0. 



(94) 
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Equation (|94p implies that we can put: 



+00 



G M (r) = e iu " T G m (u 



n— — 00 
rP 



G M K) = - / dre~ lulnT G (t) , (95) 
P Jo 



where the quantities 

_ (2n + 1)tt 



(96) 



are the fermionic Matsubara frequencies. They are such that 

rP 



dTe i{u m -u„)r = p Snm . (97) 



Since the same boundary conditions as in Eqs. (|94|) apply to Matsubara self- 
energies, we can analogously put 



E M (r) = J2 e iu " r £ M K), 



71 — — OO 

1 rP 



E (w„) = - / dre— "^E m (t). (98) 
P Jo 



This implies that 

rP 

7 M (r, t') = -i / dr!E M (r - n) • G M (n - r') 



+00 



-i/3 ^ e^ T - T ')E M K)-G M K) ; (99) 



where we have used Eq. (|97|) . Analogously. 

rP 



J M (r, r') = -i / dn G M (r - n ) • E M (n - r') 
Jo 

+00 

= -i/3 J] e iw " (T - T,) G M (w„) • E M K). (100) 



n— — 00 



Therefore, the quantities I m (t,t') = I m (t ~ t') and J M (r,r') = J M (r - t') 
depend on r and r' only via the difference (r — r'), exactly as the Matsubara 
Green functions. In the following, we will employ any of the two notations, 
according to convenience. 
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It is now useful to derive some symmetry properties of the Matsubara Green 
functions and of their self-energies. From the definition, see Eqs. ([48|) , it can be 
seen that 



which implies [see Eq.(|95 



= (-*„)]* 



(101) 



(102) 



Then, we consider the two equations of motion for the Matsubara Green func- 
tions, see among Eqs. (|87p . and compare the complex conjugate of the first one 
with the second one, with the help of Eq. fllOip and the fact that df(r— r')/dr — 
—df(r — t')/9t'. Doing this, we obtain the identity 

-'C.{T-T')]*=J^ a {T-T>l (103) 

which, written in frequency space, gives 

zZM = -K b (-" n )]*. (104) 

This implies that, in systems where = S^, £^|(t) is purely imaginary: 



E + E 

.uj n >a w„<o. 



w„>0 

= £ (e^^ fe (a; n )-e---[Sf K)]*} 
w„>0 I J 

2i Im/ge^sMK)). 



, e -^„r s 



if s™.=sr 



(105) 



8. Simplification of the action 

We now consider systematically the effective action term by term, and per- 
form some mathematical manipulations to simplify it. The general strategy 
consists in applying partial integration in order to transfer the time derivatives 
from the fields to the Green functions, at the price of receiving boundary values 
of the fields into the expression of the action, but with an important gain in the 
fact that time derivatives of the Green functions can be elaborated by means of 
the Kadanoff-Baym equations (|87|) . These are used to remove the single-particle 
matrix elements and the time derivatives of the Green functions from the effec- 
tive action, replacing them with self-energy functionals. The final goal is to get 
expressions for effective exchange parameters in terms of non-equilibrium Green 
functions and self-energies, which lend themselves naturally to approximation 
schemes. 
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8.1. Elaboration of Sqq 

We now elaborate the quantity iSqq [see the first among Eqs. (|66|) ]. In this 
case, we can directly use partial integration with respect to both integration 
variables t' and t, in order to remove the derivatives of the fields £ a 'Q(t') and 
C Q (t). We obtain 



[S QQ = E E 



dt 



te'£ Q (t) A> a ,(t,t>) + A< a ,(t,t>) Z a , Q (t') 



+>EE E 

a a' 7]—>,< 



dt< 



to 



[Cq(oo) (x£,(oo,<)G# (t,oo) - G%,(oo,t)Wf a (!,ooj) 

+ C Q (to) (Gli (to,t)W2l(t, t ) X% (t , t)Gl7 a (t, to))]ta>Q(t) 

+ (w^(t,t )Gll(t ,t) - Gf a ,(t,to)x:Uto,t)) W*o) 
-EE E («q(*o) Sg(oo)) 

a a' r/=>,< 

( Gf a , {to, t )G$ a (t , to) -Gf a , (to, oo)G# (oo, to A (^ Q (t ) ' 
V-G^(oo,t )Gl(t ,oo) G^(oo,oo)G^(oo,oo) A^'Q(°°V ' 

(106) 

where we have defined the matrix 



-^(m')-EE 



b V 



-4 



-<Wg£ + (l - ^ah) T o6 (t) 



G|i(*,t')Gfl(t',t) 



it can be proved that 



<Wi^7+7W(i') (l-K'fe') 



(107) 



(108) 



which follows from the definition, Eq. (|107l) . and Eq. ([56|) . To get an explicit 
expression for the quantity A^ a ,(t,t'), which will play an important role in our 
treatment, we first note that 



E 



1 



-<Wg7 + (l - l^ab) T ab (t) 



GfM^')Gy a (t',t) 



t/V 



= -Gf a (t', t) W^(t, if) + Xf a (t', t) Gf a , (t, f ); 



(109) 
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then, we have 
v 



J2 w ab>(t, t') tf„l, (*') = h Ri ■ x^ + e^ ■ x Al + ^ * Y ri ] / (t, t'), 



■ S At + W** ■ E^ + * S rt 



(*',*), 
(110) 

where we have defined 

xfc(t, = «(t - <5 ab - e(t' - 1) [x>°{t, o - x<* (t, f )] , 
^(t.O = *(* - t') <u + e(t - i') [w^ftO - w^(t,f )] . (in) 

Inserting Eas. (|109l) into Eq. (|107l) . and using Eqs. (|110[ ). we obtain: 

- g^, (t, t') \w& ■ + w flt • + r nt * E rt l (*', t) 

L J a'a 

. + E^ • X Ai + E^ * F^l (t, t') Gf? a {t', t). (112) 



8.2. Elaboration of Scq 

We now elaborate the quantity \Scq- The procedure is straightforward but 
a bit long, so we sketch here the main passages. We treat separately the two 
terms contributing to iScq- The first one is 



iS. 



(i) 

CQ 



dt&(*)--B (1) (*) 



(113) 



with the operator B^2 (t) defined in Eg. ([79)1 . First, we eliminate all the hopping 
terms [see Eg. ([75)1] by observing that 



Y,G< b c (t,t)T ba (t) = -X< a c (t,t)-i 



: dG< a c (t,h) 



ti=t 



£ T ab (t) G< c (t, t) = -W<° (t, t) + i 



t 1= i 



Taa' (t) G<; Q (t, t) = - (tO G<; o (tl, t) 



+ ^aa'i 



, 0G<Kfi,t) 



G$(t,t)T a , a (t) = G%{t,h) &-/ a (ti) 



<4-, 



ti=t 



ti=t 



9ti 



s 9G<i(t,tr) 



ti=t 



*i=* 



(114) 



then, we note that 

\dG<?(t,h) dG<t{tM) dG<°(t lt t) , flG<t(t ls t) 



9ti Sii 
dG<*(t,h) , dG< s (tx,t) 



0ti 



tl=i 



9ti 



0*i 



dG< s (M) = .dpS W 



dt 



dt 



(115) 



2G 



Then, the integrand contains time derivatives of both fields and £ q 'q(*)- 

It is possible to eliminate by partial integration one of such field derivatives, but 
not both. We choose to eliminate the derivative of £ a Q(t). The result is: 



+ EE [°° dt £c(t) ^al.^GtUti^ + GtUtM^-JSti) 

a a' Jt « 

- E l.Cc(oo) p S a (oo) £ qQ (oc) - CM P t(to) U)(toj\ ■ 



(t) 

(ia>Q(t) 

(116) 



The second term contributing to iScq is given by 



? (2) _ 
'CQ — 



dt 



la 



(117) 



with the operator B^,(t) defined in Eq.([5U|). Recalling the second among 
Eqs. (|66[) . and doing partial integration over dt' in order to remove the derivative 
of £ a 'Q(t'), we can write 



dCcWEE 

-J b V 



Sabig- t + ^ab) T ab {t) 



x <^ iS a , v 



00)^(00,4)^(00) - G%,(t,t)G$ a (t,t)Z a , Q (t) 



mf 



-mt 



+ 1 dt'Gf a ,{t,t')Gf a {t',t) 



Sa'b'i—+T a , b ,(t') (l-K- 



(118) 



where the index r/ E {>,<} = {+,—}. The elaboration of this term consists in 
performing partial integration, and using the sum rule 

E vGf b (t, t) Gf d {t, t) = -i5 cd G< b % t) + iS ab G<J(t, t) = 5 cdP i b (t) - 5 abP l d {t), 

(119) 

which follows from the discontinuity of the equal-times Green function, see 
Eq. (|54|) ; as a particular case, when d = a and b = c = a', this quantity equals 
— 2S aa >Pa(t). Then, one uses Eqs. (|114p , with C replaced by t or I, to remove 
the hopping terms, as well as Eq. (|112[) . and Eq. (|109p for treating the boundary 
terms. 
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After the simplifications, the total term \Scq — iS(?Q X ^CQ * s gi ven by 



/>00 pOO 



^7je(t'-t)^ a ,(t,o 



WO 



/>oc 

*EE / d CcW^aa'(*Ka' Q W 
a •'to 



^EEE^ / dtCc(*)[G^(t,oo)Xj' a (oo,t) 

a a' V *° 



^(^oo)G^(oo,i) ^ Q (oo) 



a a' V *° 

-^(*o,*)W^ B (t,to)]^g(t) 
+ E E ^c(*o) E vGf a , (to, oo) C# (oo, to) Woo) 

a a' Tj 

+ J2[^c(to)p S a (to)^aQ(to)+Cac(^)p S a(^)^Q(^)]. (120) 

a 

In Eq. (|120p we have defined 
\F aa ,{t) = S aa , [X< a s (t,t) - W< s (t,t)] 

+ i E V [X% (*, <#«(*, «) ~ ^aa' («. *) («. *)" 

= 5 aa , [X< a c (t,t) + W< c (t,t)] 

+ i { *) - *)] *) + Gl 4 -, (i, t) [w£(t, t) - W>t(t, t)] } , 

(121) 

where the passage has been done with the help of Eq. (|54|) . From Eqs. ([87j) it can 
be shown directly that 

(*> *) " X aa> (*. *) = (*> *) " (*> *) = i ^aa' (*), (122) 

while from Eq. tjiJTj) we see that 

X< c (t, t) + W< c (t, t) = - [W<°(t, t)] * + W< c (t, t) = 2i Im[w< (t, tj . 

(123) 

Therefore, we obtain the following expression: 
F aa ,{t) =2S aa ,Im\w< c (t,t)] +i^ a/ (t)G<t(t,t)+iG^(t,t)sI, a (t). (124) 
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It must be noted that, from the identity (|122p and the property (f9"Tj) . it follows 
that the Hartree-Fock self-energy is hermitean: 

poa'(*)]*=S a /„(t). (125) 

Using Eas. (|124l) . (|5t?| and (|125p . we obtain that also J- it) is hermitean: 



<F aa' (^) <F a'a(^) - 



(126) 



8.3. Elaboration of Sqc 

The term iSqc is closely related to the term iScq, which we have already 
treated. In fact, by using the properties T a b{t) = T£ a (t), as well as Eq. (|5"oT) and 



(127) 



X§aba>(t,t') = -[Xab>a>b(t',t)]* , 

it can be straightforwardly proved, starting from the definitions f|8 1 [) and (|82|) , 
that 



i-SQC 



Sqc — S, 



CQ- 



(128) 



Given Eq. (|128p , the quantity \Sqc can then be computed directly by taking the 
complex conjugate of iS C Q [see Eq. (fT2"0"]) ] with the help of Eqs.flSH]), flEE]), (HUB")) 
and CE1|. 

8.4- Elaboration of S vv 

We now elaborate the term iS vv , which is given by the sum of two terms, 
characterized by the kernel operators and D^ 2 \ defined, respectively, in 
Eq.® and Eq.®. The first term is equal to 



iS« f " drQ v (r)D^,(r)^ v (r) 

a a' J ° 

= E ^ ^ 2 ^ aa + E i P « b + ^ Pba ) 

- E (L [c m a v 03) - Cv (o)u (o)" 



£<mj(t) 



(129) 



where we have just applied partial integration to remove the derivative of £ a „ (r) . 

In order to manipulate the second term contributing to iS vv , it is convenient 
to represent the operator D^ 1 [see Eq.(|84)l] in terms of the continuous functions 
G v and G D , rather then in terms of the Matsubara Green functions G M . To 
this end, we note that 

- r>) G™1(t> - r) =6(r - r>) GZ(r - r>) G^r' r) 

+ Q(t' t) GZ(T - 7-0 G°l(r' - r). (130) 
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Wc can then write 



a a' Ja J ° 

= -EE/" d - f d ^:»EE 



0t 



To, 



rP rP 



fia'b'-g^J — T a 'b' (l — P a'b' 



EE / d -/ dr'CWEE 

a a' J ° Jt b V 

xGS(r-r')G°; o (r'-r) 



- ( i 



J* 



aft J Tab 



fia'b'-g^J ~ Ta'V ^1 — -Pa'b' 



(131) 



f2l f 2 1 *) f22^ 

After performing partial integration over dr , we can write LSi '+iSvv , 

where the first part requires one integration over the vertical branch, while the 
second part requires a double integration. The first part reads 



isg^-EEr^&ME 

a „/ JO . 



^ ab ~Q^ ~ (•"■ _ ^ab) Tib 



<C'(0) G^(0) - G^(0) G^7 a (0) f a ,„(r) 



14 / n x ^Dt 



+ GZ(T - (3) G a J a ((3 - r)t a .M - GZ(T) G v J a (-r) &,„(()) J. 

(132) 

Using Eas.([58|). we see that 

CS(0) G^(0) - G&(0) G°; Q (0) = - 5 aa ,pi, (133) 

and we further observe that 

Gg(T-f}) = -G%(T), 

GZW-T) = -G v b °(-T). (134) 



We substitute the relations in Eq. (|133l) and (|134|) into Eq. (|132p . then we per- 
form partial integration and we use the following relations, which proceed from 
Eqs.CCl]) and (£5): 



^(0)G^(0) = -(pt o ,-^)^ o ; 



(135) 
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after doing this, one needs to evaluate the quantity 



E 



1 



-<5afc7^: - (l - ~^ab) Tat 



= GZ(r' ~ r) C(r, t<) - J~V, r) (f% (r - r'), 



(136) 



in the particular case of r' = and b' = a'. The passages leading to Eq. (|136l) 
follow from the equations of motion, Eqs. (IM|) . under the condition r > t' [which 
is satisfied in the term of interest in Eq. (|132l) because r' = 0]. We have used the 
fact that, in the above expression, T ab can be safely replaced by (T ab — p,8 ab ), 
which can be immediately appreciated by looking at the first line of Eq. (|136p 
and observing that 

(l - ~f ab ) T ab f ab =(l- ^ ab ) {T ab - f i5 ab )f ab (137) 
because ^1 — P ab j S ab f a b — 0. We obtain 

f P drt av (r) + E( T <*°&a + P^ba) At |U(r)| 

... Jo b K 'Jo 

E E (T«*pla + Pta'Ta'a) f drCW L'v(t) 

EE[^w-^(°) 



i5 (2i) = _ 



Jo L 

+ E E Pta'pl'a [CM ~ CM] [ta>v(fi) " Za>v(0) 
a a' 

- E [CM pL - C(o) pL] [tavW) - U(o) 



It is now convenient to sum iSm 1 ^ with the previously computed quantity iSvv 
using the following consequence of Eqs 



(138) 



(i) 



E ( Tat >Pba - Pab T ba^j = _ * E [^^ba ( ) ~~ (0)T ba 



i J Q M a s (0, 0) - J a M a s (0, 0) + Jl [GS( n - r) + Gi s (r - n ) 



4T(0,0)-J» b (0,0) 



(139) 
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We obtain: 



LS# +i£& 1 > =i £ [/MS (0;0) _ Ja MS (0;0) j f drC(r)U(r) 

-E^[i^wi 2 -i^(°)i 

a 



x / drC(r) G^ a (-r)C(^0)-J^(0,r)Gft(r) 
Jo L 

+ E E ^a' ^'a ~ C(Q)] [^(/?) - e.',(o) 

a a' 

- E pL - C(o) pL] - Mo) 



(140) 

The second part from Eq. (JT5TJ) , requiring a double time integration, reads 

rP rP 



i^ 2) ="EE/ W dr'C(r)E 

a a' J ° J ° bb' 



Tab 



dr 



1 



I* 



-GZ(r-r')G V b ; a (r'-r) 
rP rP 



<0 



'Sa'b'g^J ~ T a >b> [ 1 - P a'b' 



^(r')e(T-T') 



EE / *7 dr'CWE 

a a' J J ° bb' 



Ot 



(ib 



Tab 



S a 'b' -q^J — T a 'b' (l — a'b' 



^(r')e(r'-T). 

(141) 



We see that, analogously to Eq. (|136p . the following holds for r > r' 

= GZ (r - r>) 7 a M l(r', r) - J^{r, T>) G V J a (r> - r 
and the following holds for r < t': 



(142) 



E G S(---')Gr; Q (r'-r) 



-Sa'b'-Q^J — Ta'b' (l — F* a i, 



GZ (r " r') Or', r) - J^(r, r') Gft^ - r) 



Ml/ 



(143) 



32 



After using these relations and integrating by parts in dr, we obtain: 



n. „> J JO h 



s n ,i - (i - p 



ab J Tab 



x {gZHt - r>) C a {r', r) - j&{T, r>) G y a (r> - r)) 



a a' 



Jo L 

C(0) GS(-r') I a M !(r', 0) + C(0) JZHO, t') GfKr')] ^ v (r'). 



(144) 

We now use the equations of motion, Eos. (jM|) . to eliminate the hopping param- 
eters. We first note that 

-$~\t) ■ J m (r, r) = [l m * E M ^] (r, r') + £ M +(r, r'), 

(r',r).V- 1 (r)= [S M ^J Mt ](r',r) + E M nr',r), (145) 



which follow from the definitions of the quantities I M and J M [see among 
Eqs. (|57|) ]. These give the following simplification: 



E 



1 



F*ab) Tab 



GZ(r-r')C a (r',r) 



Ot 

-J^(ry)G W f a (T'-T)\ 
= E ['IV, r) t- b \r) G^(t - r>) + G-l(r' - r) ^(r) J Q M >, r') 

6 

- r) K(r) GZ\ (t - r') - G^(r' - r) ^(r) J#( T> /)" 

= C(r, r>) tf a (j>, r) + J a M >, r') J^(r', r) 
+ i<5(r-r')^[/ a M ;(r,r) + J a M a Hr,r)" 

" <*£ (r - r>) * J M t] ^ {T >, t) _ G m (r - r>) E^(r', r) 

- [1^ * S^] Qa; (r, r') G^ Q (r' - r) - E£(r, r') - r), (146) 

where we have again implicitly replaced T ao with (T a b — fJ-5 ao ), which is correct as 
can be seen from the first line of Eq. (|146[) . Then, we observe that the boundary 
condition on Matsubara Green functions, Eq. (|M|) . combined with the equations 
of motion for Matsubara Green functions [among Eqs. (|87|) ], implies that 



jM (/ 3, r ) = -J b «(0,r) 

iZ(t,p) = -iZ(t,o) 



for < r < (3, 
for < r < /3. 



(147) 
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Using Eas.lfliB]). (IM1) and (fHTj) to elaborate Eq. fTIl]) , we obtain 

a a' J ° J ° 

iE /^rCMte^r) + tf(v)k m (r) 

- EE Ar [CW " 0,(0)1 [GS(-r) #*(r) - J#(-t) G&(r)j ^,(r), 

(148) 

where we have defined 

Ka* (r, r>) = C(r, r>) d(T>, r) + J a M >, r>) J^(t> , r) 

~ G£(r - r') £^( T >, r) - ^(t,t>) - r). (149) 

Finally, we take into account the fact that 7 M (r, r') = 7 M (r-r') and J M (r, r') = 
J M (r - r'). This implies that 

C S (0, 0) - Jf a S (0, 0) - ir(r, r) - J a M /(r, r) 

= 4 M a S (0) - J a M a S (0) - e(0) - j£(Q) = -J*°(0) - J Q M Q C (0) = ifila, (150) 

where, from Eas.([55| and (fTUO")) , 

1 + °° r n 

z °=2 E [CK)Ck) + CK)Erk)]=i:; (151) 

cr— b n— — CG 

the reality of this quantity is proved by using Eqs. (|102[ ) and (|104|) . 
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The total term iS vv = iSm? + iSiv + iSi 2 , 2 ^ is then written as 



i?»=-y;y dr dr'C(r) lC aa ,(T,T')+j35 aa ,5(T-T')T a ^ v (r') 

a a' J ° I 

- C(0)] [G^(-T) 7 a M l(r) - J Q M a i(-r) g£(t)] f^(r) 
+ [^,03) - ^,(0)] €» [g^(-t) 7^(r) - J Q M l(-r) (r)] | 

-E^[i^C8)i 2 -iu(o)i s 



E ^ - C,(o) p1] [UG8) - U(o) 

a 

E E pL< pL [&vW) - CM] [l>vW) - L>v(o) 

a a' 



(152) 



It is convenient to write K. in frequency space. The convolutions appearing 
in Eq. (|149[) are given by: 



+oo 



*£ M ^] (r,r') = -/3 2 E e^^'^K) • G M ^K) • S M ^K), 



n— — oo 

+ 00 



-/3 2 E e iw »( T '- T )S M ^(a; n ) ■ G M ^K) ■ J^(u n 



n— — oo 



r',r). 



(153) 



Therefore these convolutions, just as (r, r') and J^, (r, r'), depend on r and 
t' only via the difference (r — r'). We then put [cfr. Eq. (|149l) ]: 



K aa ,{T-r')= P 2 E E e^-^^^/C^K,^), 



(154) 



71 — — OO 771 — — OO 



where 

K aa >{u n ,u m ) = - [s Mi K) ■ G M ^K)1 [s M ^ m ) • G Mt (w m ) 

J aa' L 

- G M V») • S M ^K)1 fG M t( Wm ) • £ M ^ m ) 

J aa' - 

+ G^K) [s Mt K0 ■ G M n^™) ■ ^ M 



+ 



S M; K) . G M4 K) . E Mi K) G Mt (Wm) 



1 



(155) 
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From Eqs.(HnU) and (fTD3|) . it is seen that 



fcaa' {— W n , — 0J m ) — K a ' a {uJ n , W m ), 



(156) 



which in turn implies 



]Caa'(l~) = IC a ' a (T) 



(157) 



. 5. Elaboration of Sq v 
We now elaborate the term iSq v , which is given by 



/•oo pp r ~~ft , , N 

i5 Q „ = V2 e E / d * / d <Q w E E + 1 1 - p^) r «6(*) 



xG^(t,r)G 6 7 a (r,i) 
r/3 



Oa'6' 1 q Ua'b' 



Co't,(r) 



^EEE /. dT C Q (oo)G:t,(oo,r)G^(T,oo 







a a' 6' 



-CQ^o)G a i,(t ,r)G,, l a (r,t ) 

= W (l) , W (2) 



Vb'l 



<9r 



iTa'b' ( 



i - 



(158) 



where we have performed partial integration over dt. In the following, we will im- 
plicitly replace T 'b< with (T a >b>—iJ,d a >b>), analogously to what we have done in the 
elaboration of iS vv , which is correct because T a iv is multiplied by ^1 — P a 'b^j 
in Eq-dUBJ. 

Equation (| 158[) is written as the sum of two terms; the first term requires 
a double time integration (lines 1 and 2), while the second term requires one 
integration (lines 3 and 4) . The next step is the execution of partial integration 
over dr in the first term, which then requires the evaluation of the following 
quantity: 



b V 



-<W^ + (l - ^ab) T ab {t) 



Sa'b'-Q^ + Ta'b' (l — P a'b 



Gi:,{t,T)G r b ? a (T,t) 



(159) 



To simplify this object, we need the Kadanoff-Baym equations for G n and G r , 
see Eqs. (j87|) . A consequence of these equations is that Eq. ([109p is valid also 
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with the replacements /Lt — » n , n — > r , — > Y, X — > Y. Then, we get: 
GZ (t, r) = £ [GlZ> (*> i) - r) G fc ? Q (r, i)" 



= Y^,(t, r)Z r J a (r, t) + Z£,(f, r)r a r *(r, t) 

(*, r) + ■ Z*\ . (t, r) + [E^ * J^] . ( t> r) 



' =^ a , a (r, t) + [J* * a , Q (r, t) + (r, t) 



(160) 



Including the boundary contributions coming from the partial integrations, the 
first term from Eq. (|158l) is equal to 



a a' Jt l L J 



G r J a (0, t) YJ (t, 0) - Y a ,[(0, t) G a \, (t, 0) £ a , w (0) } . 



n 4 



The second term from Eq. (ll58[) is equal to 



(161) 



IS, 



(2) 



dr<< 



Cq(oo) [g£, (oo, r) Z Q ? a (r, oo) - zJJ, (oo, r) <# a (r, oo)" 



^Ee{(c q m c Q (to 

a a' L 



-G&(cx>,0)G# 8 (0,oo) G^(oo^)G^(Aoo)\ /£ o ,„(0) 



(162) 



which follows from the application of partial integration and the use of Eqs. (|87|) . 
We then sum Eqs. (|161|) and (|162|) . and we use the correspondence properties 
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between non-equilibrium Green functions listed in | Appendix A We obtain 
iS Qv =-V2iY,J2 dt dT C Q (t) OH (*, r) Uv{t) 

a a' Jt ° J ° 
a a' J *0 { 

[G&fo, t) W> a Ht, t Q ) - X<l(t Q ,t) G>J, (t, t Q )] 

[G$l(to, t) W<Ht, t ) - X>l(to, t) G< a t (t, to)] WO) | 

-^EE f^Ac Q (°o) to(c»,r)<; o (r,cx,)-^(oo ) r)^r o (r,oo) 



a a' 



■ML/ 



-fMf , 



^EE (CqM c q (*o) 



G Q < i(oo,to)G>; Q (t ,oo) G^(oo,i )G<t(i ( , x 



MP) 



(163) 



<§. 6. Elaboration of S v q 

The form of the term iS v q is analogous to that of iSq v , except for some 
substitutions. Comparing Eqs.([85|) and ((86]), we see that, in order to derive 
the expression for iS v q from that for iSq v , we must substitute in the latter: 



9. Assumption of continuous field 

By looking at the final expressions derived in the previous Section, one can 
see that the total action on the Kadanoff-Baym contour, S = Sqq + Scq + 
Sqc + S vv + Sq v + S v q can be re-written as the sum of four contributions: 

S = S cc + S hc + S ch + S b b, (164) 

where S cc includes all the contributions requiring a double integration on the 
contour, S c h and She include the contributions requiring one integration on the 
contour and depending on the values of the fields at the boundaries, while Sbb 
depends only on the values of fields at the boundaries. More precisely, we can 
parameterize the whole Kadanoff-Baym contour by means of a single parameter 
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z € [0, Zf], with Zf = 2(ioo — to) + (3 such that 

t = t + z for z G [0, ~ £o)] i G [*0j *<»] = 7-h 

t = 2t oo -t -z for z G [(too - t ), 2(too - t )] =>• t G [too,*o] =7-j 
t = z - 2(too - t ) forze [2(t oo -t ),2(t oo -t ) + /3]^re [0,/?]= 7t) , 

(165) 

and write the action as a functional of the single field £, a (z), which, with respect 
to the fields living on the various branches introduced above, is given by: 

i a {z) =e[(t 00 - 1 ) - z] e a +(t + Z) 

+ Q[z - (too - t )] e[2(too - t ) - z] e a -(2too -to-z) 

+ Q[z - 2(too - t )] £ a „(z - 2too + 2t ). (166) 

In the above equations, we have denoted as too the time value corresponding to 
the end of the branch 7+. As mentioned in Section 2, this can be taken as an 
arbitrary value, and it is usually sent to 00. In addition to the initial point z — 
and the final point z — zf, there are two other special points: z — Zf + = ^ — to 
and z = Zf- = 2(too — to), corresponding, respectively, to the end points of the 
branches 7+ and 7- . 

To proceed, we assume that the field £ a (z) is a continuous function of z on 
the whole contour, which means, in particular, that 

lim £ a (z/+ - e) - £ a (z /+ + e) = => £a+(*oo) = £ a -(*oo)i 

£->0+ L J 



lim 

£->0+ 



£«(«,_-£)-&,(«/_+£) =0 =>&,_(*«,) = £«,(<)), (167) 



so that the quantities £ a (zf+) and £ a (z/_) are uniquely defined. This assump- 
tion provides a significant simplification of the total action. It must be noted 
that, since it affects only the terms of the action labelled as She, <Scb and Sbb, 
while the magnetic interactions appear in S CCl this simplification does not affect 
the discussion of time-dependent exchange interactions, which is the focus of 
this Article, so it must not be regarded as an additional approximation. 

To exploit it, we take the total action and we separate it into the four parts 
introduced in Eq. (|164p . Then, the first among Eqs. (|167p implies £q(oo) = 0, 
which sends to zero all the terms which would require the evaluation of Green 
functions at t — >■ 00. To exploit the second among Eqs. (|167p . we transform the 
fields £c(to) and £ Q (t ) back to [f+(to) + £-(*o)]/V2 and [£+(t ) - £-(t )]/%/2, 
respectively, then we substitute £+(to) = £(0), £-(*o) = £«(0), and £„(/3) = 
£(zf). We obtain that all the boundary terms evaluated at the special points 
z — zt+ and z = Zf_ vanish, consistently with the fact that they are not 
boundary points anymore for the field £. The resulting expressions for the four 
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terms contributing to the action are: 

'A> a ,(t,t>)+A< a ,(t,ti)]e aQ (t)Z alQ (t>) 
+ [A> a , (t, t 1 ) - A< a , (t, t')] [&(* - t) C aC (t) £,a' Q {t') 

-Q(t-t')C Q (t)Za'c(t')_ 

+ iS{t- t') [Cc(t) Za'Q® + £ Q (t) &c(t)] Faa'(t)} 
poo r f3 

~EE / dr / dT 'Cv(r) K, aa ,{T-T')+(35 aa ,5{T-T')l a ^ v {t'), 

(168) 



n „' Jt JO 



^c-EE[c(o)-cm 

a a' 

■ { jf dr [GMJ (-r) #t(r) - J a M a l(-r) G^(r)] £ a ,„ ( r) 

"G<i(to,t)W>I(t,t ) - X<i(*o,*)G>^Mo)ka' Q (*) 



v^i/ 

Jt 



(169) 



i5 <* -EE ^ / d *C Q Wka > a'(i^o)G a < ; Q (to,t) - G>±(t,t )X<l(t0,t) 



Jo 



Sa>(0)-Sa>{Zf) 



(170) 



i^bb=EE[C(0) £(*/) 



a a' 



J, t 



1 1 



+ 5 a 



raa 
Paa 



t 
Paa 

-p c 

raa 



^a'(O) 



(171) 
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In the remainder of this Article, we focus on the analysis of the term S cc , given in 
Eq. (|168l) . which includes the magnetic interactions between the time-dependent 
spin fields. If one assumes that £ a (0) = £, a (zf), then the terms Sb c , Scb and 
Sbb vanish [as it can be seen by looking at Eqs. (|169l) . (|170[) and (I171[) ] . and the 
action then reduces to S cc . 



10. Exchange and twist-exchange interactions 



We now prove that the effective action for the small deviations from collinear 
spin configurations for the Hubbard model, that we have just derived, is equiv- 
alent to that of a time-dependent spin model, which means that the action can 
be written as the sum of interactions between spins. To show this, we notice 
that it is possible to transform the £ fields to unit spin vectors via the relation 



e-A ■ e-B - 1 



1 B 



lit. 



e A x e B 



(172) 



where A and B are sets of indexes, including Hubbard site-orbital, branch 
(+,—,«) and time labels. The transformation Eq. (|172j) is correct up to the 
quadratic order in the azimuthal angles 0's, which is consistent with the degree 
of accuracy of our formulation. The products of C and Q fields appearing in in 
Eq.lflggjl become: 



CqW WO = \{e*Q{t) ■ WO + i« 
£c(t)Z a >Q(t') = \{e a c{t) ■ e a , Q {t') + iu 

, e+(o-e-(o 



e a Q(t) x e a , Q (t') |, 
e a c(t) x e a , Q (t') | 



£ Q (t)L>c(t') S \{e aQ (t) ■ e a , c {t')+xu z ■ \e aQ (t) x c / (f)] } 



(173) 



where we have first undone the Keldysh transform from (£q,£c) to (C+i^-) m 
order to apply Eq. (|172|) . and then we have defined 

e aC (t)= [e a+ (t) + e a -{t)]/V2, 

e aQ (t) = [e 0+ (t) - ea-(t)] /\/2. (174) 
After this transformation, we can write S cc as the sum of two contributions: 

*5cc = ^exchange ~t~ ^twist ; (t75) 
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where 

I poo poo ( 

iSoxchangc = lEE / d< / d *' 1 2ie ac(0 ' e a 'Q(t') 
4 a a' -^O ^to [ 

• ~S(t - i') Re \T aa , (<)] + 0(f - i) Im [,4> , (i, - .A< , (t, 
+ e aQ (t) ■ e a , Q (t')Re [A> a , (t, t') + A< a , (t, t')} j 



a a 7 



to 



dr' 



e ot) (r) • e a / 1 ,(r') - 1 /C aa /(r - r') 



(176) 



is the exchange part (including all the terms depending on scalar products &a 
es), and 

I r^o poo ( 

iStwtot = ~ A u z ■ EE / d W d M 2ie «c(t) x e a , Q (t') 

4 a a' "'to "/to { 

S(t - Im [Jw (t)] + 6(t' - t)Re [U>, (i, i') - .A< , (i, 0] 

3«g(t) x e a , Q (f')Im [.A>,(t,t') +A< a ,(t,t')] j 
1 



+ • E E y t d< y o dT e «Q (*) x e °'-( T ) (*« r ) - g li> (*> r ) 

-7 iu *'EE / dr f P to J ea*(T)xe a , v (T J )K aa ,(T-T') (177) 

4 a a' J ° J ° 

is an additional part, including all the terms depending on vector products 
x efl), that corresponds to a new kind of magnetic interaction, which we 
call twist exchange. In Eq. (|176|) we have already removed an obviously null 

term equal to e av (r) ■ e a ' t ,(r') — 1 8 aa i8{r — r')/3I a , and we have exploited the 
symmetry of e a ±(t) ■ e a >±(t') under the interchange of both a a' and t <-> t' , 
combined with the symmetry properties of the quantities A^ ,(t,t') and J- aa '(t) 
under the same transformations. Similarly, in Eq. (|177|) we have removed an 
obviously null term equal to e av (r) x e a i v {r') S aa /S(T — r')/3I Q , and we have 
exploited the antisymmetry of e a ±(t) x e a t±{t') under the interchange of both 
a «-> a' and t -R> t' . 

From the transformation Eq. (|172[) . it is seen that also an action term Sg* 
appears, depending on the squares of the azimuthal angles A . This anisotropy 
term, however, is identically zero as a consequence of the symmetry properties 
of the kernel matrices. We give the detailed proof of this fact in | Appendix B| 
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Therefore, we correctly find that the system does not exhibit any anisotropy, 
consistently with the fact that we are not including spin-orbit, and the action 
is entirely written in terms of interactions between spins. 

The term Stwist * s particularly interesting, and deserves special attention. 
While at a first sight it seems to represent an effective Dzyaloshinskii-Moriya 
interaction, we know that this cannot be the right interpretation, since our 
system does not include spin-orbit. On the other hand, the action of a non- 
relativistic system must be invariant under rotation of all spins of the same angle 
about an arbitrary axis. Therefore, the allowed effective interactions between 
spins must be of the forms oc a\ ■ er 2 , or oc (<ti x er 2 ) • er 3 , or products of these 
combinations. While the first form is commonly identified as the exchange 
interaction, the second one is less known. To the best of our knowledge, the 



first discussion of a similar interaction was given by Bogoliubov [46( in 1949, in 
the context of the theory of the polar model [47| (a precursor of the Hubbard 
model). His treatment is based on what we would call today a half- filled single- 
band Hubbard model (in equilibrium) , and a three-spin interaction of the form 
oc (<Xi x <x 2 ) • er 3 appears at the third order in the expansion of the Hamiltonian 
in powers of T/U. However, in that case the corresponding term vanished by 
symmetry Here, we are not making a perturbation expansion in terms of T/U, 
since in our expansion the quantities which are assumed to be small are the 
angles of deviation of the spins from the equilibrium configuration. Therefore, 
we receive all terms allowed by symmetries, irrespectively of their dependence 
on the Hubbard parameters. Since our expansion is up to the second order in the 
deviation angles, it corresponds to interactions oc Seri ■ <5<r 2 and, without loss of 
generality, oc (<5<xi x <5<t 2 ) • (T3. In the second form, only cr\ and cr 2 are rotated, 
while the third spin vector is pointing along the original quantization axis, which 
is u z in our model. Therefore, since the third vector is fixed, the interaction has 
the form oc (ei x e 2 ) -u Zl which is precisely what we find. Therefore, despite its 
resemblance with the Dzyaloshinskii-Moriya interaction, we conclude that the 
contribution S'twist describes a three-spin non-relativistic interaction; the unit 
vector u z is to be interpreted as the quantization axis of the third spin, which 
is not rotated in the expansion. Note that the Dzyaloshinskii-Moriya vector is 
not parallel to u z [l8[ . Differently from the case studied by Bogoliubov, we are 
here treating a multiband Hubbard model in non-equilibrium, and in our case 
this interaction does not vanish in general. 



11. Equilibrium 

Before treating the general non-equilibrium case, which is the main goal of 
this Article, we here specialize the discussion to the equilibrium regime. This 
had already been studied, e.g., in Ref. ljl, 17 1 with different methods, and we 
now show that we are able to reproduce the known results. 

If the system is in equilibrium, the general treatment pursued so far is 
overkill. In fact, consider the expression of the time-dependent expectation value 
of an operator O, as given in Eq.(|5]). In the equilibrium limit, i.e., H{t) = Hq, 
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one has U(t, t') = cxp[— i(t — t')Ho]. Exploiting the cyclic property of the trace, 
as well as [N, Hq] = 0, one obtains the equilibrium expression: 



Tr 



0(t) 



Tr 



-P(h -hn) 



Tr 


U v (t - ij3,t )6 


Tr 


' e -/3(H -vN)~ 





(178) 



which shows that the equilibrium expectation value of any operator can be 
evaluated by just considering the vertical branch of the Kadanoff-Baym contour. 
In the path-integral formulation, this results in an effective partition function 
given by 



(179) 



where «S E Q = —iS vv . Only iS vv is therefore needed to compute the equilibrium 
properties, and in this Section we study this portion of the action in order to 
recover the known equilibrium results. 

11.1. Equilibrium action 

The equilibrium action (excluding the boundary terms) is given by: 

S EQ [e(r)] 4 EE T dr f &T>K aa ,{T-T'){e av (T)-e alv (T>)-l 

4 a a' J ° J ° 1 

+ iu z ■ \e av { T ) x e a , v {r') }. (180) 

This action is a functional of the r-dependent helds e av (r). For the particular 
choice of constant fields, i.e., 



Lv(t) = U Vtg[0,/3], 



(181) 



which is equivalent to neglecting the so-called quantum fluctuations (see, e.g., 
the discussion in Chapter 6 of Ref.[48j), after integrating over dT and dr' the 
quantity JC a a'{T — t') we obtain the action 



S EQ [e] =^ 4 EE E K «'K. W ») e«»'ea'.-l + m 2 '(e„xe a /„ 



a a/ n— — oo 



From Eq. (|156p , one can show that 



+oo 



+oo 



^ ICaa'(u 



7i •> ) 



(182) 



(183) 



71 — — OO 
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and by using this relation together with the symmetry of (e av ■ e a / v — lj and 
the antisymmetry of e av x e a > v under a O a' , we obtain 



av v 



1) Jaa> + (e av x e /„) ■ ii 2 y aa ' , (184) 



where we have introduced the equilibrium exchange parameters 

+00 



-Re 

4 



2J ^Caa'(w„,( 



(185) 



and the equilibrium twist- exchange parameters 

+00 



-ilm 

4 



2J ?Coo'(w„,£J n ) 



n— — 00 



(186) 



Evidently, from Eq. (|183p it follows that J aa > = Ja'a and y aa > = ~y a 'a- 

We can make a further step by re- writing Eqs. (ll85[) and (|186l) in the r- 
representation, which is achieved by inverting the Fourier transforms: 



+00 ^ ,p ,p ,p „/3 +00 

V] }C aa '(iJ n ,U} n ) = — / dri / dr 2 / dr 3 / dr 4 V] e 

P Jo Jo Jo Jo „__™ 



-iw„(ri+r2+r 3 +r4) 



r 2 ; 



r 3 



T4 



X M ^(r 2 ) G M ^(r 3 ) £ Mt (74) 

- aa' - 

+ GS(n) [s Mt (r 2 ) • G M nr 3 ) • E M t(r 4 )' 



/? 4 Jo 



+ OO 



dn / dr 2 ^ 



-icj„(ri+r 2 ) 



GS(n)E5S(7a) + ES(n)GK(n) 



Ml, 



n— — oc 



Now, the sums of the form 

+00 



+00 



J2 e-^ T = cos(u; n r) 



(187) 



(188) 



n— — 00 



n— — 00 



are real. Moreover, since the system that we are studying is non-relativistic 
(NR), we can choose a basis set of real single-particle wave functions, which 
implies that T a f, is real, hence the driving operator on the LHS of the Matsubara 
equations [the last two among Eqs. (|87|l ] is real, as well as the equilibrium pair 
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correlation function p ab . From Eqs. (|55|) . it then follows that the Matsubara 
Green function G^J(t) is purely imaginary close to the discontinuity point t = 0, 
and will therefore stay imaginary on the whole domain —fi < r < j3. This 
also implies that the self-energy £^(r) be imaginary. With this choice of the 
basis set, it must be noted that from Eq.flTUI]) it follows that Gf b (r) = Gf a (r), 
and from Eq. (fTU5)) it follows that £^,(t) = E*(t). Therefore, the quantity 
appearing in Eq. (|187l) is real, and from Eq. (1186[) we have y aa > = 0, consistently 
with the discussion by Bogoliubov 

11.2. Local single-band case in equilibrium 

To establish the correspondence with the results of Ref. , we now study 
the equilibrium action (|184[) in the local single-band (LSB) approximation. Here 
local means that the self-energy is assumed to be diagonal in site index, and 
because of the fact that each site carries a single orbital (single band), we put 

Ej£>„) ^ foE^K). (189) 
Under this assumption, from f|155p we obtain 

+oo +00 
n=— oo n= — oo 

--f £ [G^K)sr t K) + srK)G^K)_, 

n— — oo 

(190) 

where E^ s = (E^ - S^) /2. The term proportional to 5 00 / appearing in 
Eq. (|190p docs not contribute to the action, as it can be seen by looking at 
Eq. (|184l) . hence it can be immediately suppressed. Eliminating this term, in 
the LSB case we can write the effective exchange parameters just as 

J aa , L I B - £ EfK) G$K) S«K) G^K), (191) 

n— — oo 

which is precisely the result obtained in Ref.|16l| with a different method. Our 
analysis provides the more general expression Eq. (|185|) . which can be applied 
in the case of a non-local (non-diagonal in site index) self-energy, as well as 
the interpretation of the results of Ref.[l6| in the more general framework of 
the path-integral formulation. Our action includes also terms depending on 
r-dependent spin vectors (quantum fluctuations), which contribute to the com- 
plete treatment of the equilibrium properties. 

12. Non-equilibrium 

We now consider the full non-equilibrium action, given by Eqs. (ll76p and 
P77|) . We notice that both the exchange and the twist-exchange contributions 
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are written as sums of several different terms: f ) a term S RR requiring a double 
integration on the real-time branch, which is related to the real-time evolution 
of the system, 2) a term S RI requiring one integration on the real-time branch 
and one integration on the imaginary-time branch, which describes the influence 
of the initial correlations on the evolution of the system, and 3) a term S 11 
requiring a double integration on the imaginary-time branch, which describes 
the system in equilibrium. 

In turn, the first term can be written as the sum of three different contribu- 
tions. The term involving products of two Q spin vectors is the source of noise 



32l . l37j . which has contributions of both exchange and twist-exchange forms. 



The pertinent parameters are: 

J?Zr(t,t>) EE IRC [A> a ,(t,t')+A< a ,(t,t')] , 

yrJ sc (t,t>) ee -Jim [A> a ,(t,t') +A< a ,(t,t')} . (192) 

The term involving the products of a C and a Q spin vectors can be in- 
terpreted by first going back to the representation in terms of branch vec- 
tors, e + and e_. One has 2e aC (t) ■ e a , Q (t') = e a+ (t) ■ e a > + (t') - e _(t) • 

e a'-(i') + e a _(£) ■ e a i + (t') — e a+ (t) ■ e a /_(i') , where the first (second) pair 

of square brackets enclose a term which is symmetric (antisymmetric) with re- 
spect to the tranformation (a,t) ■<-» (a',t'). Analogously, 2e a c(t) x e a >Q(t') = 

e a+ (t)xe a > + (t')-e a -(t)xe a ,_(t') + e a _(i) x e a <+(t') - e a+ {t) x e a /_(i') , 

where the term enclosed by the first (second) pair of square brackets is anti- 
symmetric (symmetric) under (a,t) ■(-> (a',t'). Then, exploiting the symme- 
tries provided by Eqs. (|I08p and (|I26[) . we can write the resulting coefficient of 



e a -\-(t)-e a i+(t')— e a -(t)-e a '-(t') in Eq. (|I76l) . which we call the non- equilibrium 
exchange parameter, as 

J aa ,(t,t') = -~{5(t - t')MFaa'(t)} ~ Im[^l> ,(t,t') -A< a ,(t,t')] j, 

(193) 



and wc can write the coefficient of the term 



e a+ (t) x e a , + (if ) - e _ (t) x e a ,-(t') 
in Eq. (|177|) . which we call the non- equilibrium twist- exchange parameter, as 

y aa , (t, ee I {s(t - 1') im \T aa , (t)\ + Slgn(f 2 ~ f/) ite [A> a , (t, t>) - A< a , (t, t')]\. 

(194) 

It can be immediately seen that J a a'(t,t') — J a 'a(t\t), while y aa '(t,t') — 
—ya'a(t',t), analogously to the properties of the equilibrium parameters J aa i 
and y aa '- The structures of the non-equilibrium parameters, as defined in 
Eqs. (|I93[) and (|I94I) . are analogous to those of the equilibrium quantities, but 
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the Matsubara Green functions appearing in the latter are replaced here by 
real-time non-equilibrium Green functions. The last contribution to S RR is rep- 
resented by additional terms which are related to the causality structure of the 
action: we define 

AJ aa ,(t,t') = ~lm[A> al (t,t') - A< a ,(t,t')] , (195) 

and 

Ay aa , (t, t') = - iRc [A> a , (t, t') - A<, (t, t')] . (196) 

These two parameters have anomalous symmetry properties: AJ aa i (i, t') — 
— AJ a ' a (t\ t) and Ay aa i(t,t') = Ay a ' a (t',t). Hence, they can be called, respec- 
tively, anomalous exchange and anomalous twist-exchange parameters. In terms 
of the quantities introduced above, the exchange and twist-exchange parts of 
S can be written as 

^change =EE/ dt dt'\e aQ (t)-e a , Q (ti)J™r(t,t') 
-2ie aC {t)-e a , Q (t') 

and 

/>oo /"OO [ 

i^twfst =EE / dt / dt '{ X e °'^')] • U ^a°J SC (t,t') 

a a' Jt Jt { 

- 2i [e aC (t) x e a ' Q (t')} ■ u z 

It is apparent that the computation of the parameters entering Eqs. (|197|) 
and (|198p , as well as of the remaining ones from Eqs. (|176l) and (|177l) . would be 
cumbersome in the most general case. In order to obtain a tractable problem, as 
well as to gain insight into the physical processes contributing to non-equilibrium 
magnetic interactions, we shall discuss two approximations. The first one is the 
time- dependent Hartree-Fock (HF) approximation, which consists in neglecting 
memory effects in the electronic self-energies. The second one is the assump- 
tion of slow spin dynamics (SSD), i.e., that the dynamics of the spin degrees 
of freedom be much slower than the dynamics of the electronic degrees of free- 
dom. The two approximations are independent because they involve different 
quantities (the spin vectors in the former case, the self-energies in the latter 
case), therefore it is possible to apply either just one of them, or both. In the 
following, we discuss them separately. 

12.1. Time- dependent Hartree-Fock approximation 

We now consider, for the general non-equilibrium case, the time-dependent 
Hartree-Fock (HF) approximation [3l[. This amounts to neglecting all self- 
energies contributions having time arguments corresponding to different points 



J aa ,(t,t r )+AJ aa ,(t,t') 



(197) 



y aa '{t,t') + Ay aa ,{t,t') 



(198) 
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on the Kadanoff-Baym contour. In other words, from Egs. ([571) we neglect S > , 
E < , £ n and £ r , and we keep only the equal-time singular parts of the self- 
energies, i.e., from Eqs. dS^l) we obtain 

t') = F T, A (t, t') = F S(t - t') E(t), (199) 
and, analogously, the Matsubara self-energy becomes 

E M (r) = F /35(r)S M , (200) 

where we have introduced the factor (3 just for convenience. This approximation 
allows to compute trivially all the integral convolutions between Green functions 
and self-energies. Physically, it corresponds to neglecting memory effects in the 
self-energies. The technicalities implied by this approximation are discussed in 
|Appcndix D[ while here we will present the results for the main quantities of 
interest. In the case in which the self-energy can be assumed to be diagonal in 
Hubbard indexes, the action can be further simplified as shown in | Appendix F[ 
The expression for the exchange and twist-exchange parameters can be sim- 
plified by means of the representation of Green functions in terms of time- 
dependent occupation numbers and of the spectral function, which is defined 
as 

A(t,t') = i [G R (t,t') - G A (t,t')] = i [G>(t,t') - G < (t,t')] . (201) 

In equilibrium conditions the Hamiltonian is time-independent, therefore the 
Green functions G>(t,t') depend only on (t — t'), hence it is possible to take 
their Fourier transform with respect to the relative time, 

/+oo 
dte^Gf a m (t). (202) 
-OO 

It can be shown (see the proof in [Appendix E[ ) that the following relation holds: 

G< EQ (,),i/(,)CH, (203) 

where 

f(u) = wrp. 7 (204) 

is the Fermi distribution, and A^®(lj) is the Fourier transform of the equilib- 
rium spectral function. 

In non-equilibrium the Green functions Gf a (t,t'), as well as the spectral 
function, can be written as functions of t r = t — t' and T = (t + t')/2. Taking 
the Fourier transform with respect to t r leads to the following definitions: 

G<(cj,T) = if{cj,T)A ba (uj,T), (205) 
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where 



+00 



G<>,T) = / d(t-t')e^-^G< a (t,t% 



(206) 



and f(w,T) and Ab a (uj,T) are, respectively, T-dependent occupation number 
and spectral function. Equation (|205l) is effectively the definition of f(ui,T). 
The definition of the spectral function, Eq. (|201|) . then implies G^ a (uJ,T) = 
i[f(u,T)-l]A ba (u>,T). 

The transformation from the Green functions G> to the spectral function 
and the time-dependent occupation number is particularly useful in the HF 
regime. Using Eqs. (|D.2l) and (ID.3|) . we obtain the following expressions: 



A< a ,(t,t')TA> a ,(t,t>) 



HF 



dw 
2^ 



dui' 

1 

2tt 



j4 



T 



A^{lu',T) 



A±(u>, T) • Wt - %)] \a^(oj', T) ■ E f (t + % 



T 



■ A 1 {lo,T)-T, 1 [T 



where t r =t-t' and T = (t + t')/2, and 



f-(u,T) = f(u,T)-f(w',T), 

f+(u, T) = 2/(w, T) /(«', T) - /(w, T) - JV, T); 



(207) 



(208) 



*W (*) 4'>> *) + 4a' (w, t) SL(<)1 (209) 



where we have removed the term proportional to 5 aa / coming from T aa i [see 
Eq. (|D.3[) ]. because it gives a null contribution to the action. In fact, the 
quantity T a a'{t) enters the action in the product with S(t — t') and either 



(*)-e«/+(f) 



-(*) 



-(*') 



e a+ (t) x e a > + (t') - e„_(i) x e a /_(t') 



Both combinations of spin vectors vanish identically for a = a' and t = t' , 
therefore the components of T a a' proportional to S aa > can be suppressed. 

The quantities of interest can be written down using Eqs. (|207[) and (|209[) . 
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For example the exchange parameters, from Eq. (|193p . become 



HF 



J aa > (t, t') n = -5{t r ) Re<J / — f(co, t) p; o , (t) AJ,>, t) 



+ Sign{U) Im< 



+ °° dw 

-oo 27r 
t r 



doj' 

( 

2tt 



-i{uj—uj')t r 



E T (T-|)-AV,r) 



f(u>,T)-f(cj',T) 



+ 



sMt + ^).^( w ,t).eMt-^ 



7=T-l 



(210) 



S7ow spin dynamics 

To implement the condition of slow spin dynamics (SSD), we look back to 
the term iS RR of the action, which is given by Eqs. flWj) and (fHgl) . While the 
spin vector e a q is associated with the fast quantum fluctuations of the local 
magnetization, we can recognize the classical magnetization in the spin vector 
e a c, see e.g. Ref . [39i] . It is then possible to assume that the ec vectors be 
slowly- varying functions of time with respect to the rapid electronic processes, 
which are embedded in the exchange and twist-exchange parameters. Therefore, 
one can expand 



g(t,t')e aC (t)ttg(t,t') 



(t - 1) 2 

e a c(t') + (f - t')e aC (t') + y - ' e aC (t') 



(211) 



where g(t, t') is an electronic-related function (such as the exchange parameters). 
No approximation can be done, instead, on the quantum vectors eQ, since their 
dependence on time is fluctuatory and not under control. When only the term 
of zero order in (t — t') is kept in the expansion (|211[) . one obtains 



2i EE 



dt 



la 



dt'e aC (t) ■ e a , Q {t') J aa >(t,t') + AJ aa ,(t,t') 



2i / dt ' e *c(t') ■ e alQ {t') \jaa>(t')+AJ aa ,{t') 



(212) 
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and 



dt / dt' [e aC (t) x e a , Q {t')\ ■ u z \y aa ,(t,t') +Ay aa , (t,t') 

to J to 



2i EE 

a a' 

' -2i I dt ' [ e "c(t') x e a , Q {t')] ■ u z y aa ,(t') + Ay aa ,(t') 



(213) 



where 



POO POO 

Jaa>{t')= / dtj aa ,(t,t'), AJ aa ,{t') = / dtAJ aa ,{t,t'), 
/•oo />oo 

y aa .{t')= dty aa ,(t,t'), Ay aa ,(t')= dtAy aa ,(t,t') 

J to J to 



(214) 



are the parameters that account for spin precession due to magnetic interac- 
tions within the Hubbard system. The terms in the expansion (|212l) depending 
linearly and quadratically on (t — t') generate, respectively, spin damping and 
spin nutation |39j . 

The SSD approximation may be useful for computational purposes because, 
if the expansion (|211l) is truncated at any order in (t — t'), it gives an effective 
action where the classical spin vectors e a c (and their derivatives) are evaluated 
at the same time t' as the quantum spin vectors e a <g. The price to pay is that 
the parameters of the magnetic interactions, depending on two times (t,t'), 
must be integrated over dt on the whole real axis. A technical simplification 
can be achieved by using the sum rule Eq. (IB. 14[) , combined with the symmetry 
properties (|108|) and (|126[) : for example, the exchange parameters can be written 
as: 



Jaa'{t') 



\ |Re[.7w(t')] - \ J™ dUign{t-i?)lm[A> a ,{tJ) --A^CM')] } 



1 



1 







dtIm[A> a ,(t,t')-A< a ,(t,t') 



(215) 



For large values of t' , this representation avoids the need of computing an in- 
tegral over dt from to to t' , replacing it with an integral over the imaginary 
branch. 

It is then clear that the further level of simplification consists in just ap- 
plying both HF and SSD. This can be done by computing Eqs. (|214p with the 
approximate expressions for A, T and Q given in |Appendix D| [Eqs. (|D.2|) . (|D.3j) . 



13. Dynamical spin stiffness 

The spectrum of spin waves [43j | in ferromagnets, in the non-relativistic case 
and in the absence of external fields, can be computed directly from the spatial 
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Fourier transform J„ of J aa i , as 



u> q = -^{j -J q ), (216) 



where M is the magnetic moment per ion. In the long-wavelength limit, the 
dispersion uj q is quadratic, 



q->0 



a p 



where a, /3 € {x,y,z} and the tensor -D Q( 3 is known as the spin stiffness. This 
tensor is a relevant quantity for comparison with experiments, since it can be 
determined from measurements of t he sp in- wave spectrum. Furthermore, for the 
equilibrium case it has been proved [49| that, if the self-energy is local, then the 
expression for the spin stiffness computed by applying Eq. (|191j) and Eq. (|216jl is 
exact, despite the approximation of neglecting the vertex in the Dyson equation 
for the two-particle Green functions, Eq. (l4"Tl) . 

However, our analysis shows that when an external field is applied, the sit- 
uation is more complicated. First, there is a non-zero time-dependent twist- 
exchange interaction, which may in principle alter the spin-wave spectrum. 
Moreover, even in regimes where this effect is negligible, the presence of a time- 
dependent external field leads to a modification of the concept of spin-wave 
itself. Only in the SSD case, corresponding to adiabatic spin dynamics, we can 
think of generalizing Eq. (|216l) as 

« -[ Jo (t)-J q (t)\, (218) 

where J q (t) is the spatial Fourier transform of J a a'{t) given by Eqs. (|214p . Equa- 
tion (|218|) assumes that the magnetic moment per ion M remains constant and 
that J aa ' (t) depends spatially only on the lattice vector R aa > = R a — R a > , as 
in equilibrium. This is justified if the external time-dependent field can be as- 
sumed to be almost uniform over a portion of the sample large with respect to 
the interatomic distances, a condition which is usually met in experiments on 
laser-induced spin dynamics, where the typical diameter of uniformity of a laser 
beam is of the order of some tens of /im. 

We shall therefore consider for simplicity the case of a single-band Hub- 
bard model, and assume that the non-equilibrium twist-exchange interactions 
be negligible. Under these conditions, wc apply the HF approximation, obtain- 



ing effectively the HFD approximation discussed in Appendix F To determine 
the spin stiffness, we put 



^(M') = ^£e ife ^'G7 (t,t'), (219) 



G 



where n is the number of k vectors in the first Brillouin zone (1BZ). We use 
Eq. (|F.3[) for the exchange parameters, imposing also the SSD approximation, 



53 



and we directly get 

Jaa>{t') = -Y J j* R ™'J q (t'), (220) 



n 



where 

(221) 

In the small \q\ limit, we can expand 

(222) 



but the term linear in q a does not contribute to Eq. (|221[) because Gk = G_fc. 
Computing Eq. (|218[) with the help of Eq. (|22ip . we find indeed that the disper- 
sion is quadratic, with a spin stiffness tensor given by 

/>oo 

A*^— 2M&E / dt'sign(t'-t)Z S m S (t') 
n a Jt 

1 dGT(t>,t) dGf(t,t>) 
Details on this derivation are given in |Appcndix G| 



14. Conclusions 



In this Article we have provided for the first time a method to compute the 
parameters characterizing the evolution of the low-energy magnetic interactions 
in a general realistic system (a multi-band Hubbard model) under a generic 
time-dependent electrostatic perturbation. While confining the analysis to the 
low-energy sector (small spin deviations from the equilibrium configuration) is 
crucial to obtain a tractable problem, there is no restriction on the time-scale 
of the perturbing field, which allows us to describe ultrafast excitations. Im- 
portantly, our method requires only the computation of Green functions which 
are diagonal in spin, which is expected to be a major simplification in view of 
numerical implementation. 

Our findings show that the evolution of the system is a highly non-trivial pro- 
cess, involving complicated interactions between the macroscopic magnetization 
and quantum fluctuations. In particular we find the emergence of a new mag- 
netic interaction, the twist exchange, that is formally similar to a Dzyaloshinskii- 
Moriya interaction, despite the absence of spin-orbit in the system that we have 
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treated. This term can be interpreted, generalizing Bogoliubov's work [46(, as 
a three-spin interaction. Moreover, we have derived expressions for noise con- 
tributions arising from the interaction between quantum fluctuations. Finally, 
we have discussed relevant approximations (time-dependent Hartree-Fock and 
slow spin dynamics), and we have showed how to compute the time-dependent 
spin stiffness tensor, which may be relevant for comparison with experiments. 

Future developments may involve the incorporation of spin-orbit interaction 
and external magnetic fields, and the derivation of effective equations of mo- 
tion for the magnetization. The numerical implementation of the formalism 
described here, which is best done in the framework of time-dependent dynam- 
ical mean- field theory E3-52|, combined with employing appropriate hopping 
and interaction parameters, should allow to simulate realistic systems under 
ultrafast (femtosecond) perturbations. 
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Appendix A. Kadanoff-Baym functions with particular arguments 

Some relations hold between different Kadanoff-Baym Green functions, for 
certain values of their arguments. In this Appendix we list such relations. We 
first observe that 



(A.l) 



Gl, a (t,0) =i(^ a ,(t)) = G< a (t,t ), 



and 



G r aa , (0,t) = -i{ % (0) ij>l, (f ) } = -i ( e (t) 



Gl, a {t, p) = i = i (Vv(<) = -G>, a (t, t ), 



G r ba (r,t ) = -HMr)ti 



GZ(r) 



(A.2) 



Gl(t ,r) =i(#(T)&) = [&(0)#(t)]) = GZ(-r). (A.3) 



55 



As particular cases, 

G^^.to) = -i(&09)$) = = -G^(to,*o) = -ipta, 

G^to, /3) - i (#00) 4) = i = -G> (to, to) = -ipta + i<W (A.4) 

As a consequence of Eq. (IA.1[) , 

n (t,o), 

(0,*). 
(A.5) 



r n (t,o) = 


£ fl -G n + £ n *G M " 


(t,0) = 


£ R • G<~ 


(t,*o) + 


E n *G M " 


r r (o,t) = 


G r • £ A + G M * S r " 


(0,t) = 


G> • 


(t ,*) + 


G M *£ r " 



As a consequence of Eq. (|A.2|) and Eq. 



^ (*,/3) = 



F r (/3,t) = 



£ fl -G n + £ n *G R 



£ R -G> 



(*,*o) 



G • £ + G * £ 



G< -£ y 



(*o,*)- 



M) = - 
£ n *G 

GM) = " 
G n *E 



£« G> 



(Mo)- 



£ ★ G 



(t,t ) = -W^t.io), 



G < • 



(*o,f) 



G * £ 



(t ,t) = -X < (t ,t). 



(f,0) 

(0,t) 
(A.6) 



As a consequence of Eq. (jA.3|l , 



W<(t,to) 



X>(t ,t) = 



E^ • G< + E n * G r + £< ■ G A 



(Mo) 



: R G< 



{t,t ) 



(t,o) = r n (t,o), 



G> • E A + G n ★ E r + G 



(to,*) 



G ; 



(to,*) 



G * £ 



(0,t)=F r (0,t), (A.7) 



where we have used Gf a (t',to) — as well as G R a (t ,t') — because f' > to- 
Then, it also holds that 



Z n (t ,r) = 
Z r (r,t ) = 



G fi -£ n + G n *£ M 
£ r -G A + £ M *G r 



(*o,r) 
(r,t ) ■ 



G * £ 

E M * G w 



(0,r) = J M (-r), 
(r,0) = / M (r). (A.8) 
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Appendix B. Proof of Sg2 = 

The action term produced by the transformation Eq. (|172l) and depending 
on the squares of the azimuthal angles is equal to: 



+ /V[U<,(t,i') ->4> ,(* J tO + ^.(*'»*) - ^.(f.*) 



i / dr 



^(t,r)+^(t,r) 



t4- 



«E/ ^Wi^+E/ dT / [^/(r-r / )+^ a (r / -T) 



(B.l) 

We shall now prove that Sg2 = 0. 

We start by computing the last line of Eq. (|B.l|) . which contains terms related 
to the vertical branch of the Kadanoff-Baym contour. We have: 

fP r 1 + °° r 

/ dr' /C oq /(t - t') +/C / (t' - t) =/3 3 V] /C OQ / (u n , uj n ) + K, a , a (uj n , u 7 

J ° 1 J n=-oo L 

(B.2) 

where we have used Eqs. (|154p and (p?Tl) . Then, we need to sum the RHS of 
Eq. (|B.2[) over a': we use (| 155[) as well as a sum rule valid for the Matsubara 
Green functions, which reads 

/3 2 G M ^K) ■ S MS K) • G Mff K) = -G MS (w„); (B.3) 

we present the full derivation of this sum rule in |Appcndix C| One then obtains: 



a' 

and 

a' 

which imply 



1 r 

T 2 



S>„)'G^k) + G Mt K)-S Mt k) , (B.4) 



l r 

T 2 



G^( Wn ).E M H Wll ) + ^K).G«t( Wn ) , (B.5) 



^ ^aa'(Wn,W„) + /C a ' a (w n , CJ„) = —2pX a , 



(B.6) 



as can be seen by comparison with Eq. (|151[) . Therefore, the last line in Eq. (|B.l 
is equal to zero. 
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To proceed, we need to elaborate the following sum and integral: 



E / df A< a ,(t,t')-A> a ,(t,t')+A> a (t',t)-A< a (t',t) 



(B.7) 



To this aim, it is convenient to write the quantities A>, from Eq. (|107|) . as: 



(Bi 



and 



^a(*'.*)=E 



6' L 



• [ - W^t(*, Gf t) + CfJ, (t, X&(f , t)] . (B.9) 

In this way, when one computes the sum over a' in Eq. (|B.7I) , the terms depending 
on the hopping parameters give zero by symmetry, and the remaining terms are 
derivatives with respect to t', which can then be integrated trivially. With the 
help of Eq. (|122|) and the definition (|124[) . we then obtain 



E / d*' A< a ,(t,t')-A> a ,(t,t')+A> a (t',t)-A<, a (t',t) 

= - J E yaa'(t) +Fa>a(t) +E ^ (t, t )G<l(t Q , t) - G>£ (t, t )X<l(t , t) 



W<J (t, h)G>l(to, t) + G<Z(t, t )X>°(t , t) 



We now evaluate the following sum and integral: 



(B.10) 



-iE /^^,(*,T) + ^,(t,T)l =-iEE f^GZ{t,T). (B.ll) 



a' u 



To do this, it is convenient to write the quantities Q^,, from Eq. (|159|) . as 
Gil, (*, t) = - E [ - Y ™> (*> r) G£(t, t) + Gil, (*, (r, *)" 



(B.12) 
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and then proceed analogously to the previous computation of Eq. (|B.7[) . We 
obtain 

-*EE [ ^QZ(t,r) 

= -i E E [ - Y «*'^ G ^a{P, t) + Gil, (t, P) Y7a(0, t) 



a' a 



+ Y%(t, 0) G r J a (0, t) - Gil, (*, 0) Y2(0, t) 
E E [ " W aJ(t, to) G<l(t , t) + G>° (t, to) X<l(t , t) 



a' er 



+ W<J(t,t ) G>l(t ,t) - G<:,(t,t )X>l(t ,t) , (B.13) 



where for the last passage we have used Eqs. dA.ll) . (IA.7I) . (IA.2|) . (|A.6|) . Collect- 
ing all the terms, one obtains 



/3 



i / dr 
'o 



dt' 



n 



A< a ,(t,t')-A> a ,(t,t') + A> a (t',t)-A< a (t',t) 



= 0, 



(B.14) 



which, together with Ea. (IB.6l) . proves that Sg2 = 0. 



Appendix C. Sum rule for Matsubara Green functions 



The proof of Sgi = 0, presented in |Appcndix B[ requires a sum rule that 
we will now derive. We start from writing the KB equations of motion for the 
Matsubara Green functions, among Eqs. (|87|) . in the frequency domain. To do 
this, we multiply both sides of the KB equations by e _1 "' l ' T_r ) and we apply 



l P S T d(r - t'), where < n < 0. After performing the integration with the 
help of Eqs. flU}, (HSJ) , ^ and (fTUD|) . we obtain 

1 



1 



G M >„) • g~\ Un ) = /3G Mcr (ui n ) ■ E M<T (w„), 



where we have put 



iffafc 1 ^") = -i^ndab ~ {Tab ~ fJ-5 a b)- 



(C.l) 



(C.2) 



We now multiply both sides of the second among Eqs. (|C.ip by G Mcr (w n ) from 
the right, obtaining 

G MCT K) • g-\w n ) ■ G Mff K) = - Q G M "{Lo n ) - /3G Mff K) • £ Mct K) ■ G M <>„), 



(C.3) 
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then we use the first among Eqs. (|C.l|) to eliminate g 1 (w n ) in the LHS of 
Eq. (|C.3l) . and we obtain the required sum rule: 



^-/3E Mff K)-G Mff K) 



G (uj„) 

G Ma (uj n ) — G Ms (u> n ) 



-g ~ fiG {uj n ) -E (w n ) 



G M5 K) 



- f5G Ma {ijjn) ■ E m<t (cj„) — E Mff (cj„) • G M9 (io n ) 



/3 2 G M ^K) • E MS K) • G M <>„) = -G MS ( W „). 



(C.4) 



Appendix D. Technical consequences of the Hartree-Fock approxi- 
mation 

We here list the identities which follow from adopting the Hartree-Fock (HF) 
approximation. From Eqs. (j87|) . using Eqs. (|199l) and (|200[) we obtain (matrix 
notation): 



HF 



(t,t') = ^(t)-G^(t,t% 



HF r^fiVr „"l„1 /, s HF ^<T 



(t ) t , ) = G^(t,0'S'(t') l 



Y n<T (i,r) = [E fl ' CT • G^] (t,r) = E (*) • G n<7 (t,r), 
Y r<T (r,i) = F [G rCT -E^] (r,i) = F G rCT (r,i) -E"(t), 
Z^,r) ^ F [G^*S Mct ] (t,r) == -iG^(t,r) • E~ M ", 
Z rCT (r,t) ?= F [E MCT ^G rCT ] (r,t) = F -iE Ma • G r *(r,t). 



(D.l) 



Inserting these approximate identities into Eqs. (|112l) . (|124j) and (|160[) . we ob- 
tain: 



G^(i,i') 'sV)-G $t (* , ,*)-S T (t) 



G^(i, f') • eV)1 fG^(t', t) • E T (t) 

J aa' - 



(D.2) 



HF 1 



E T (<) • G <1 ~(t, i) + E ; (i) ■ G<^, t) + G <t (t, i) ■ E f (i) 



+ G<+(t, t) • E (t) + i EL, (t) G<t(t, t) + i G<j; (t, t) Ej, a (t), 

J aa 

(D.3) 
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GZ(t,r)^i \>f(l).G"(t. \ GJjr.n 
+ iG^(t,r) 



S Mff .G rff (r,t).r(i) 



G r5 (r^).E CT (i) 



(D.4) 



Appendix E. Occupation numbers representation in equilibrium 

For the sake of completeness, we here prove Eq. (|203p . In equilibrium, H (t) = 
Hq, and G > (ti,t 2 ) becomes 

Gl EQ (t) = {e-^-^^^e-^l} , (E.l) 

where Z = Tr |exp -/3 (jl - fiNj j , and t = t 1 -t 2 . We put K = H - fiN, 
and we observe that 



hence we can write 

Gf a EQ (t) = -^Tr {e^te^e-^} ; 
analogously, one obtains 

{e- a *$e it *&e-''*} . 



G 6 < EQ (*) = 11 



2 



-Tr 



(E.2) 
(E.3) 

(E.4) 



In both expressions, we explicitly express the trace by employing a complete set 
of many-body eigenstates \n) of K (with eigenvalues K n ), and we also insert 
the spectral decomposition of the identity between -0J and ipb- We obtain 



-itfj. 



C(*) = -VEE e " 



i(iC„-if m ) p -^if„ 



G 6 < EQ (*) 



-ii/i 



(K n -K m )-PK n 



The Fourier transforms with respect to t are 

27ri 



Gta Q H = -r E E^ w + K n K ™ + M) e" 



■01 



06 



n), (E.5) 
ra\ . (E.6) 



x ( n 



■0i 



(E.7) 



and Gjf Q EQ H = -exp [/3 (w - /x)] Gjf Q EQ (w). Equations ((2031) and flU]) follow 
immediately. ■ 
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Appendix F. Hartree-Fock approximation with diagonal self-energy 



In the Hartree-Fock approximation, if it is possible to further assume that 
the self-energy be diagonal in Hubbard site-orbital indexes (which we denote as 
HFD approximation), the expressions that we have derived for the dynamical 
parameters of the magnetic interactions undergo a remarkable simplification. 
This approximation is somehow analogous to the LSB approximation that we 
have discussed in the equilibrium case and, analogously, it may be relevant 
only in the special case of a single-band Hubbard model (with a local Hartree- 
Fock self-energy). From the identities following from the HF approximation, see 
Appendix D| with the additional assumption of diagonal self-energy, we get: 



xf:,(t,t') n i u Gf:,(t,t')r a ,(t% 

^(r,t) H = D G^(r,t)S:,(t), 

Z aa>(t, T ) = -lGw(*,T)i; a , , 
HFD . =Mer 



= -iCC(M); (F.l) 



it follows that 



H =° -^(t')G^(t\t)W a (t)G^t% 
T aa ,(t) -AiS aa ^ S a (t)G<!(t,t) = 46 aa ,X S a (t) p s a (t), 

QZ^t) H l D 4i^(t)G^(t,r)E^! S G r J a (r,t). (F.2) 

Since T aa < (t) is now proportional to S aa i , it does not contribute to the exchange 
and twist-exchange actions, hence it can be removed from the expressions of 
the parameters. After eliminating the irrelevant term, the expression for the 
exchange parameters becomes very simple: 

J aal (t,t') H ^sign(t-t')El,(t')^ (F.3) 

rj a 

It must be noted that the part multiplying sign(i — t') is identically zero for 
t — i', so one has J^ a a'(t,t) H = D 0, with no ambiguities related to the definition 
of the function sign(x) . The twist-exchange parameters in HFD conditions are 
given by 

y aa > (t, t>) H ^ D -sign(t - t')E S a , (t')Z S a (t)ij>J> G%, (t, t>) Gl? a (f , t). 

(F.4) 



G2 



In the occupation numbers representation, 



. ,f\ HFD . S/ t r \-S/ t r \ 1 f +M d, f M d W ' 

J aa ,(t,t) = - slg n(t r )X a (T+-)Z a ,[T--) -J_^ -J_^ — 



x e 



-i(uj—uj')t r 



f(u>, T) - /(a/, T)\ Ka> (w, T) A* T), 

(F.5) 



and 



x e 



(F.6) 

where t r =t-t' and T = (t + i')/ 2 - 



Appendix G. Derivation of the spin stiffness 



We here add some details on the derivation of Eq. (|223p . Using Eq. (|22ip . we 
directly obtain 



poo 

D a p(t)=—J2vYl / di'sign(i'-i)£ S (i)xV 



' X! rj,, ' T "' 



dk a dkp 



(G.l) 



To remove the second derivatives of the Green functions, we first put 



Q-qa ° ^k 



d 



dk a dkp dk a 



d 



<ji tt; — <j\ 

: /3 



rfc 9fc«" fc 



<9fc Q 9fc,g 



(G.2) 



then we transform the sum over the k vectors in Eq. ljG.lj) into the corresponding 
integral over the volume V of the 1BZ, and we use 



d 

v dk a 



d 

fc dk B fc 



dV = 



in ■ 



i) 



&k dk~0 k 



dS a , 



(G.3) 



where dS a is the projection along the a direction of the oriented infinitesimal 
element of surface of the 1BZ at k. Separating the contributions with k a > 
and k a < 0, we see that they cancel with each other, because Gk = G-k while 
dS a changes sign under k a — > —k a . Therefore, the integral of the first term of 
the RHS of Eq. (|G.2[) is zero. The remaining term gives Eq. (|223|) . 
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